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Final  Technical  Report 


Contract  AFOSR  F4962096-1  -0476 
9/30/96-9/30/98 

“WDM  Laser  Sources  for  the  Defense  University  (Testbed)  Research 

Internet  Program  (DUTRIP)” 

I.  Research 

The  objective  of  this  project  was  to  define  specifications  of  a  hybrid  integrated 
wavelength  division  multiplexed  laser  source,  to  identify  the  appropriate  vendor,  and  to 
negotiate  an  acceptable  price  structure. 

II.  Design  Specifications 

In  consultation  with  a  number  of  user  groups  and  the  DUTRIP  program  at  University  of 
Maryland  (PI:  Prof.  Mario  Dagenais)  we  have  determined  the  following  set  of 
specifications  for  the  four  wavelength  WDM  array: 

1)  Emission  wavelengths  of  DFB  lasers: 

channel  1:  1549.32  nm 

channel  2:  1552.52  nm  (the  reference  of  193.1  THz.) 

channel  3:  1555.75  nm 

channel  4:  1558.99  nm 

All  channel  wavelengths  to  be  accurate  to  ±  0.3  nm* 

2)  SMSR  >  30  dB  under  40  mA  peak-to-peak  modulation  and  8.2  dB  extinction 
ratio  (SONET  OC-48  spec.) 

3)  Threshold  current  .  <  30  mA 

4)  External  efficiency  .  >0.2  mW/mA 

5)  Fundamental  transverse  mode  operation  up  to  Ij^c  =100  mA 

6)  Power  coupled  into  .  >  +6.0  dBm 

single  mode  fiber  @100  mA 

7)  Modulation  bandwidth  .  2.5  Gb/s  ** 


8)  Four  ECL  inputs  to  drivers:  ECL,  50  Q 

9)  Four  single  mode  outputs,  optical  isolator  in  each  laser  package. 

10)  Back  facet  monitor  in  each  laser  package 

1 1)  Front  panel  setting  of  laser  bias  current  and  temperature  for  each  laser. 

12)  Front  panel  indicator  lights  to  indicate  operation  of  each  laser 

*  All  wavelengths  and  spectral  properties  measured  at  a  chip  power  output  of  5  mW.  The 
wavelength  may  be  trimmed  with  a  TC  cooler,  as  long  as  other  specs  are  maintained 
**  modulation  bandwidth  is  limited  by  the  driver  chip. 

Additional  Considerations 

The  vendor  will  (a)  provide  available  test  and  reliability  data,  and  (b)  establish  the  device 
code  and  make  it  available  to  other  customers  at  a  fixed  price 

III.  Results 

These  specifications  were  discussed  with  a  number  of  vendors.  Ortel  Corporation  was 
identified  as  a  low  cost  supplier  capable  of  satisfying  all  of  the  above  specifications.  Ortel 
Corporation  fabricated  and  delivered  52  arrays  of  WDM  lasers,  each  operating  at  four 
wavelengths,  as  specified  above.  We  have  established  a  testing  methodology  for  these 
WDM  arrays  and  conducted  extensive  test.  The  results  are  summarized  in  an  attached 
paper,  presented  at  the  1999  SPIE  Photonics  West  in  San  Jose. 

We  have  investigated  wavelength  accuracy  and  stability  of  WDM  transmitters 
delivered  to  us  by  Ortel  Corporation.  Transmitter  performance  was  tested  under  DC  bias 
and  modulation  at  data  rates  as  high  as  2.5Gbps.  The  wavelength  accuracy  of  hybrid 
devices  was  found  to  be  considerably  better  than  ±0.1  nm.  A  considerable  fraction  of 
devices  tested  showed  small  wavelength,  less  than  0.02  nm,  and  power,  less  than  0.1  dB, 
excursions  due  to  digitization  noise  of  control  circuitry.  Transmission  experiments 
through  100  km  of  standard  single  mode  fiber  show  low  chirp  penalties  and,  in  case  of 
monolithic  arrays,  the  absence  of  cross-talk.  Long  term  transmission  experiments,  up  to 
1000  hrs,  through  wavelength-selective  elements  simulating  WDM  systems  do  not  show 


any  penalties  due  to  the  wavelength  or  power  instabilities.  All  the  devices  tested  in  our 
experiments  are  suitable  for  WDM  systems  research.  This  work  has  allowed  us  to  < 
conclude  that  the  transmitters  delivered  by  Ortel  met  or  exceeded  our  specifications. 

In  order  to  facilitate  testing  we  have  constructed  a  number  of  test-beds,  some  of 
which  included  passive  devices  such  as  multiplexers  and  demultiplexers.  These  passive 
devices  had  to  be  independently  characterized.  In  the  process  we  realized  that  the  current 
modeling  of  passive  devices  could  be  considerably  improved. 

The  design  of  waveguide  components  and  optical  circuits,  requires  detailed 
knowledge  of  light-guiding  parameters,  field  distributions,  transfer  -efficiencies,  etc.,  in 
complex  device  structures  and  devices  coupled  into  systems.  Much  of  the  previous  work 
done  in  this  field  concentrated  on  two  aspects  of  this  problem.  First,  several  powerful 
methods  have  been  developed  to  evaluate  eigenvalues  and  eigenfunctions  for  specific 
waveguide  structures.  These  are  the  finite-difference  or  multi-grain  method,  the  finite- 
element  or  variational  method,  the  effective  index  methods,  the  boundary  element  or 
integral-equation  method  and  the  specific  matrix  implementation  of  it  known  as  the 
Galerkin  method,  the  WKB  method,  numerous  numerical  methods,  etc. 

Second,  after  the  eigenvalues  and  eigenfunctions  are  calculated,  their  propagation 
through  an  optical  circuit  must  be  accounted  for.  The  propagation  aspects  makes  the 
calculation  of  field  distributions  a  complicated  and  time-consuming  task.  Much  of  the 
past  effort  has  been  devoted  to  the  computation  of  field  distributions  and  propagation 
characteristics  of  waveguide  devices  in  which  optical  energy  couples  between  the  input 
and  output  ports.  The  best  known  examples  are  the  work  of  Dragone,  and  Henry  and 
coworkers,  who  used  finite  Fourier  transform  propagation  method.  Very  little  attention 
has  been  given  in  the  past  to  coupling  between  waveguides  and  free  space  optics. 

The  previous  work  usually  considered  field  distributions  in  optical  waveguide 
structures  in  terms  of  scalar  electrodynamics.  However,  the  pioneering  work  on  radar, 


antennas  and  microwave  waveguides,  has  demonstrated  the  general  usefulness  of 
describing  electromagnetic  fields  in  terms  of  vector  electrodynamics.  The  significance  of 
this  comes  from  the  fact  that  the  optical  waveguide  systems  are  described  not  only  by  the 
wave  equation  for  the  fields,  but  also  by  the  boundary  conditions  which  are  vectorial  in 
nature.  Moreover,  in  optoelectronic  devices  the  distribution  of  polarization,  a  vector-like 
quantity,  must  also  be  considered.  Consequently,  a  vectorial  formulation  might  be 
expected  to  have  an  advantage  in  computation  and  design  of  optical  waveguide  systems 
and  in  modeling  of  the  interface  between  waveguides  and  free  space. 

We  have  developed  a  new  transfer-amplitude  method  to  describe  field 
distributions  and  efficiencies  of  optical  waveguide  systems  based  on  vector 
electrodynamics.  The  starting  point  of  the  analysis  is  Huygens'  picture  of  propagation  of 
time-harmonic  electromagnetic  fields  in  complex  waveguide  structures  as  expressed  by 
Love's  formulation  (  A.  E.  H.  Love,  Phil  Trans.  Royal  Soc.  London  A,  v.197,  pp.1-45, 
1901)  of  the  vectorial  Green's  theorem.  This  approach  allows  us  to  divide  a  complex 
optical  waveguide  structure  into  a  set  of  primitive  blocks  and  calculate  a  transfer 
amplitude  and  a  transfer  function  for  each  block  separately.  We  demonstrate  that  Love's 
expressions  for  equivalent  surface  current  densities  lead  to  correct  results  for  the  field 
distributions  in  primitive  blocks.  It  is  then  possible  to  calculate  partial  and  total  transfer 
functions  and  corresponding  device  efficiency  with  high  accuracy. 

The  representation  of  a  complex  waveguide  structure  by  a  set  of  principal  surfaces 
separating  primitive  blocks  with  known  equivalent  surface  currents  and  propagation 
functions  was  shown  to  be  a  very  useful  method.  This  method  provides  a  possibility  of 
constructing  the  total  system  response  as  a  convolution  of  responses  of  its  primitive 
blocks.  It  is  a  very  close  analogue  of  the  method  of  geometrical  optics  with  its  principal 
surfaces  and  primitive  blocks.  The  distinction  is  that  an  exact  propagation  function, 
including  diffraction  effects,  is  used  here  to  describe  each  primitive  block,  instead  of  the 
simple  ray  approximation.  These  ideas  were  already  applied  to  the  evaluation  of 
Dragone's  router,  a  waveguide  problem  of  considerable  complexity.  These  ideas  can  be 


also  applied  to  modeling  of  3D  problems  which  include  both  waveguide  and  free  space 
propagation,  eg.  coupling  of  lasers,  with  an  arbitrary  far  field  distribution,  to  fibers 
through  lenses  and  lens  systems.  The  model  will  include  arbitrary  source  polarization  and 
will  allow  for  calculating  the  sensitivity  of  coupling  efficiency  to  component 
misalignment. 

Two  papers  describing  these  ideas  in  detail  are  attached  to  this  report. 
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Abstract 

This  paper  reviews  performance  of  hybrid  and  monolithically  integrated  WDM 
transmitter  arrays  based  on  directly  modulated  2. 5 Gbps  lasers,  with  a  focus  on  the 
wavelength  accuracy  and  stability  under  normal  operating  conditions.  We  also  consider 
power  stability,  chromatic  dispersion  penalties,  and  the  channel  cross-talk.  Large 
numbers  of  four-wavelength  devices  were  obtained  and  evaluated  under  a  program 
designed  to  provide  university-based  system  researchers  with  advanced  WDM 
components  We  show  that  multi-wavelength  laser  arrays  capable  of  high-performance 
out-of-the-box  operation  can  be  now  produced  for  research-type  WDM  systems. 


Introduction 

Multiple- wavelength  laser  sources  are  the  fundamental  component  of  WDM  systems  and 
considerable  research  and  development  effort  has  been  devoted  to  their  design  and 
fabrication.  The  research  effort,  in  particular,  has  focused  on  monolithically  integrated 
multiple-wavelength  laser  arrays,  widely  believed  to  be  a  superior  source,  particularly  in 
applications  involving  many  wavelengths.  However,  while  practically  all  of  current 
commercial  WDM  systems  rely  on  discrete  lasers  with  pre-selected  wavelengths, 
properties  of  such  hybrid-integrated  sources  remain  poorly  documented.  This  divergence 
of  research  and  industrial  practice  underscores  the  difficulty  of  producing  stable, 
monolithic,  multiple-wavelength  arrays  capable  of  meeting  WDM  system  requirements. 
However,  with  WDM  systems  moving  towards  higher  channel  density  and  higher  bit 
rates  the  virtues  of  the  hybrid  vs  the  monolithic  approach  need  to  be  reassessed. 

In  the  design  of  multiple  wavelength  laser  arrays,  either  hybrid  or  monolithic,  one  has  to 
consider  a  range  of  issues: 

•  Gain  material:  a  typical  WDM  system  has  a  20  to  30  nm  operating  range  in  the  1 .55 
micron  region  which  is  defined  by  the  Er-doped  fiber  amplifier  (EDFA)  gain  window. 
Lasers  based  on  quantum  well  structures  that  exhibit  continuously  tunable  gain  region 
of  more  than  200  nm  might  be  the  structures  of  choice  [1].  More  recently,  there  has 


been  also  increasing  interest  in  the  use  of  fiber  lasers  in  WDM  applications,  as  their 
gain  matches  that  of  the  EDFA  [2]. 

•  Wavelength  tuning/locking:  High  quality  laser  sources  can  be  readily  produced  using 
a  distributed  feedback  grating  (DFB)  or  a  distributed  Bragg  reflector  structure  [3], 
Currently  most  semiconductor  laser  arrays  utilize  DFB  structures.  With  increasing 
number  of  WDM  channels,  the  channel  spacing  becomes  smaller  and  the  wavelength 
registration  accuracy  requirements  become  more  stringent.  Current  dense-WDM 
systems  call  for  channel  spacing  of  50  GHz,  with  a  GHz  wavelength  accuracy.  In 
order  to  operate  such  laser  arrays  under  different  environmental  and  operating 
conditions,  careful  control  of  the  temperature  [4]  and  power  are  required.  A  number  of 
laser  designs  incorporating  additional  build-in  electrodes  for  wavelength  tuning  have 
been  proposed  [5].  In  a  sampled  grating  configuration  [5a],  more  than  60  nm  tuning 
was  demonstrated.  Wavelength  stabilization  through  external  frequency-selective 
feedback,  provided  through  an  external  bulk  grating,  fiber  grating,  or  a  WDM 
demultiplexer,  is  also  being  investigated  [6].  In  addition,  active  frequency  locking 
using  solid  state  etalons,  or  other  types  of  wavelength  lockers,  is  under  study  for  the 
long  term  wavelength  stabilization  of  single  and  multiple-frequency  laser  sources 

[6a]’  .  .  „ 

•  Power  combining:  The  use  of  simple  3  dB  couplers  to  combine  to  combine  2  laser 

sources  results  in  the  power  loss,  for  each  laser,  of  N  x  3  dB !  With  the  use  of  a  star 
coupler  the  loss  is  still  10  x  logN.  In  a  40  channel  WDM  system,  this  amounts  to  a 
loss  of  16  dB!  As  the  channel  number  increases,  it  becomes  desirable  to  use 
wavelength  multiplexers  to  combine  the  outputs  of  a  multiple  wavelength  source.  The 
typical  interference  filter  or  fiber  grating  WDM  multiplexer  for  8  to  16  channels  has  a 
fiber  to  fiber  loss  of  2-3  dB.  As  the  channel  count  increases,  both  the  insertion  loss 
and  the  cost  of  multiplexers  increase  rapidly.  For  large  channel  count  applications 
integrated  WDM  multiplexers,  such  as  arrayed  waveguide  grating  (AWG),  become  a 
viable  choice  [7].  The  16-channel  AWG  typically  have  a  fiber-to-fiber  insertion  loss 
of  4-6  dB  and  the  loss  increases  only  slightly  for  larger  channel  counts.  However, 
impressive  research  results  to  the  contrary,  such  devices  cannot  be  readily 
incorporated  into  a  monolithic  array. 

The  optical  feedback  is  another  serious  issue  related  to  power  combining.  While 
dealing  with  discrete  lasers,  one  can  use  isolators  to  eliminate  the  optical  feedback. 
However,  isolators  suitable  for  integrated  laser  arrays,  in  particular  high  channel- 
count  laser  arrays,  simply  do  not  exist  and  the  optical  feedback  remains  a  very 
difficult  problem. 

•  Data  modulation :  Current  WDM  systems  operate  at  2.5  Gbps  and  the  data  rate  is 
moving  toward  10  Gbps.  To  provide  high  quality  optical  signal  at  2.5  Gbps  and  higher 
rates,  the  data  should  be  encoded  through  external  modulation.  For  individual  lasers, 
and  hybrid  WDM  sources,  this  can  be  readily  accomplished  via  external  EO 
modulators.  To  modulate  integrated  laser  arrays,  a  monolithically  integrated  on-chip 
modulator  is  desirable.  There  has  been  considerable  progress  on  integrating 
electroabsorption  modulators  with  DFB  lasers  [9,  9a].  Electronic  laser  driver  arrays 
have  also  been  demonstrated  at  data  rate  of  2.5  Gbps  and  beyond  [10,10a].  In  hybrid 
integration  of  electronic  driver  chips  to  a  laser  array,  via  wire  bonding,  electrical 
cross-talk  can  be  an  issue. 
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The  early  integrated  multiple-wavelength  arrays  were  based  on  DFB  and  DBR  two- 
section  lasers  allowing  for  electrical  tuning  [11-12].  Large  laser  arrays  were  made  based 
on  bulk  [13]  and  strained  quantum  well  [14-16]  active  layer  structures.  The  strained 
structures  produce  lower  threshold  lasers  with  narrow  linewidth  and  lower  thermal 
crosstalk.  These  devices  required  either  high  resolution  lithography  to  produce  the  proper 
DFB  grating  pitch  or  needed  repeated  holographic  exposures  [12].  Introduction  of  the 
phase  mask  was  a  major  advancement  in  the  lithographic  mask  definition  [19].  E-beam 
lithography  was  also  used  to  prepare  highly  accurate  gratings  with  fine  pitch  shifts.  To 
further  adjust  the  wavelength,  individual  on-chip  heaters  were  introduced  [17-18]. 

Early  WDM  experiments  combined  the  output  of  laser  arrays  into  a  single  fiber  with 
the  use  of  elaborate  bulk  multiplexers  [20].  Miniature  diffraction  gratings  and  hybrid 
micro-optics  were  introduced  to  simplify  this  process  [21].  An  on-chip  power  combining 
element  for  a  three-laser  array  was  first  reported  by  Koren  et  al.  [12].  An  integrated 
optical  amplifier  was  also  provided  to  compensate  for  the  combiner  loss.  An  NxM  star 
coupler  was  later  used  to  combine  the  output  power  of  a  20- wavelength  laser  array  [22]. 
In  a  separate  effort,  electroabsorption  modulators  were  inserted  between  the  lasers  and 
the  star  coupler  to  permit  high  data-rate  modulation  of  individual  lasers  [23].  In  both 
cases,  an  optical  amplifier  at  the  output  was  used  to  compensate  for  the  power  loss. 

Although  WDM  multiplexer  is  a  more  efficient  power  combiner  for  a  large  channel- 
count  laser  array,  it  is  quite  difficult  to  match  the  MUX  passband  to  the  laser  wavelength 
in  a  monolithic  chip.  However,  the  MUX  can  be  used  as  a  wavelength  locker  to  generate 
multiple  wavelength  outputs.  The  MAGIC  (multistripe  array  grating  integrated  cavity) 
laser  [24]  was  the  first  demonstration  of  this  concept  utilizing  a  curved  mirror  grating  as 
the  wavelength  selective  element.  The  MAGIC  laser  not  only  produces  multiple 
wavelength  output,  it  also  combines  them  into  a  single  output  channel  [25].  The 
wavelength  locking  idea  was  also  implemented  with  the  arrayed  waveguide  grating 
(AWG)  WDM  multiplexer  [26-29].  The  multi- wavelength  router  laser  can  be  also  used  as 
a  tunable  laser  source.  A  digitally  tunable  laser  source  with  build-in  electroabsorption 
modulator  has  also  been  reported  [31-32]. 

Finally  we  note  the  wavelength  selectable  fiber  ring  lasers  using  Er-doped  fiber 
amplifier  as  the  gain  medium  and  an  AWG  for  wavelength  selection  [2,  33].  Since  the 
laser  cavity  is  long,  there  are  typically  103  to  104  axial  modes  within  the  3  dB  passband  of 
the  AWG.  The  laser  can  potentially  operate  multimode  and  thus  introduce  excess 
intermodal  optical  beat  noise.  A  recent  paper  reported  the  use  of  a  semiconductor  Fabry- 
Perot  optical  amplifier  as  an  intracavity  narrow  band  filter  to  stabilize  laser  oscillation  in 
a  single  axial  mode  [34]. 


Experimental  Results 

The  hybrid  integrated  laser  sources  used  in  this  work  were  provided  by  the  Ortel 
Corp  (transmitter  model  10348A).  The  design  goal  was  to  produce  two  four-wavelength 
sources,  each  on  an  ITU  grid  of  3.2  nm  (400  GHz),  with  the  shortest  reference 
wavelength  of  1549.32  nm  and  1550.92  nm  respectively  for  each  WDM  source.  Two 
such  sources  can  be  easily  combined  into  an  eight  wavelength  source  or  a  source 
operating  on  a  200  GHz  grid  with  wavelengths  varying  from  1549.32  to  1560.61  nm. 
Each  transmitter  module  consisted  of  four  independent  plug-in  units,  each  containing  an 
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ECL  compatible  laser  driver,  a  distributed  feedback  laser,  an  optical  isolator,  and  an 
FC/PC  fiber  optic  connector.  Individual  lasers  were  independently  temperature  and 
power  controlled.  LED  indicators  were  provided  to  indicate  the  status  of  control  circuitry. 
At  the  operating  point,  the  fiber-coupled  power  output  from  each  lasers  was  set  to  be 
greater  than  2  mW.  Four  separate  plug-ins  were  inserted  into  a  rack  mountable  chassis  to 
provide  a  DC  power  supply  for  the  necessary  bias  voltages  to  each  transmitter.  Over  one 
hundred  such  transmitters  were  delivered  and  tested. 

The  monolithically  integrated  laser  arrays  were  obtained  from  Nortel  Corp.  Each 
source  consisted  of  a  monolithic  four- wavelength  laser  array  fiber-pigtailed  to  an  array  of 
4  single-mode  fiber.  The  light  from  the  4-laser  array  was  simultaneously  imaged  through 
an  optical  system  composed  of  3  successive  lenses  and  a  single  optical  isolator.  An  ECL 
compatible  laser  driver  was  used  to  drive  the  laser  at  2.5  Gbps.  The  laser  array  was 
mounted  on  a  single  temperature  controller,  with  a  single  control  circuit  for  the  four  laser 
array. .  The  fiber-pigtailed  power  in  each  fiber  is  more  than  1  mW.  The  wavelength  of  the 
laser  array  varied  between  1552.52  to  1557.36  nm  in  steps  of  200  GHz.  The  DC  bias 
current  for  each  laser  could  be  controlled  individually. 


Wavelength  (nm) 


Fig.  1  A  set  often  spectra  obtained  from  hybrid  four-laser  WDM  sources  provided  by  the 
Ortel  Corporation. 

WDM  laser  transmitter  modules  were  tested  for  wavelength  and  power  stability.  For  the 
initial  DC  testing  of  hybrid  arrays,  light  output  from  individual  lasers  was  combined 
using  a  commercial  1x4  combiner,  displayed  on  an  optical  spectrum  analyzer  with  a 
resolution  of  0.8  A,  and  stored  in  a  computer.  A  spectrum  from  each  array  was  obtained 
after  a  five  minutes  warm-up  time,  to  allow  the  power  and  temperature  control  circuits  to 
stabilize.  All  the  arrays  were  tested  directly  after  shipment  from  the  manufacturer. 
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without  any  adjustments  to  control  circuits.  A  set  of  spectra  obtained  on  ten  hybrid 
sources  is  shown  shown  in  Fig.l,  on  a  logarithmic  scale.  The  side  mode  suppression 
better  than  38  dB  was  measured  for  all  the  channels.  The  output  power  of  the  lasers,  as 
preset  by  the  manufacturer,  showed  variation  of  a  few  tenths  of  a  dBm.  The  spacing 
between  the  channels  was  nominally  400  GHz,  as  specified  by  the  manufacturer,  and  all 
the  hybrid  WDM  units  tested  passed  the  manufacturing  wavelength  accuracy 
specification  of  ±0.1  nm.  The  control  wavelength  of  each  laser  can  be  front  panel 
adjusted  to  an  accuracy  of  order  0.01  nm. 

Test  results  obtained  on  12  monolithically  integrated  arrays  are  shown  in  Table 
I.  No  temperature  adjustment  was  performed  on  the  array.  By  slightly  adjusting  the 
temperature  of  the  arrays,  a  wavelength  accuracy  of  ±  0.15  nm  with  respect  to  the  ITU 
grid  can  be  obtained. 

Table  I 


Serial 


Wavelength  error  versus  corresponding 
ITU  standard  wavelength  (nm) 


Number 

ITU  1552.52 

nm 

ITU  1554.13 

nm 

ITU  1555.75 

nm 

ITU  1557.36 

nm 

980803-1 

-0.15 

0.03 

iheehhh 

-0.01 

980803-2 

-0.14 

0.00 

0.00 

980803-3 

-0.15 

0.13 

0.02 

0.03 

980810-1 

0.13 

0.11 

0.04 

980810-2 

0.07 

-0.12 

0.03 

980810-3 

-0.17 

0.13 

-0.15 

0.03 

0.08 

-0.08 

980826-2 

-0.12 

0.04 

-0.12 

■ 

980826-3 

-0.20 

0.05 

-0.07 

0.02 

-0.10 

-0.28 

-0.12 

0.10 

-0.13 

-0.09 

980817-2 

-0.16 

0.11 

0.06 

-0.05 

More  detailed  temporal  wavelength  stability  measurements,  for  lasers  under 
modulation  at  a  rate  of  2.5  Gbps,  are  shown  in  Fig.  2.  These  measurements  were  carried 
out  on  a  monolithically  integrated  array  using  a  wavelength  meter  with  an  absolute 
wavelength  accuracy  of  ±  .005  nm  and  a  display  resolution  of  0.001  nm.  In  each  panel 
the  ITU  reference  wavelength  is  indicated  by  a  continuous  line,  the  data  is  shown  as 
dotted  lines.  Three  of  the  lasers  shown  here  are  stable  with  time,  to  within  0.015  nm.  The 
fourth  device  exhibits  small  frequency  jumps,  on  the  order  of  0.02  nm,  representing  the 
digitization  noise  of  the  temperature  control  circuit.  Similar  noise  events  are  observed  in 
hybrid-integrated  lasers. 

Fig.  3  illustrates  the  temporal  power  stability  of  monolithically  integrated 
devices.  As  expected,  there  is  a  good  correlation  between  the  wavelength  and  power 
jumps.  In  lasers  which  do  not  exhibit  digitization  noise  events,  the  power  level  is 
maintained  constant  to  within  ~0.1  dB.  However,  each  noise  event  results  in  a  power 
jump  of  about  0.25  dB.  Similar  jumps  occur  in  hybrid-integrated  arrays 
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S/N  N9808 17-202 


Figure  2.  Wavelength  stability  of  a  Nortel  monolithic  laser  array 


Figure  3:  Temporal  power  stability  of  a  Nortel  monolithic  laser  array 


Bit- error-rate  measurements  through  100  km  of  standard  fiber  were  carried  out 
on  both  types  of  arrays  in  order  to  evaluate  the  chirp  penalty.  This  is  an  important 
consideration  in  closely-spaced  WDM  systems.  Even  though  DFB  lasers  are  known  to 
have  a  narrow  linewidth  under  CW  conditions,  these  lasers  will  exhibit  frequency 
chirping  under  modulation.  This  frequency  chirping  will  lead  to  a  dynamic  broadening  of 
the  laser  linewidth  and,  because  of  chromatic  dispersion  in  a  conventional  fiber  (i.e. 
Coming  single  mode  SMF-28  fiber),  will  lead  to  a  penalty  in  receiver  sensitivity,  as 
measured  in  a  bit-error-rate  (BER)  measurement.  It  is  therefore  important  to  characterize 
WDM  transmitters  in  a  transmission  measurement.  BER  measurements  were  performed 
on  both 


Fig.  4  Schematic  diagram  of  a  bit-error-rate  measurement  performed  on  both  hybrid  and 
monolithic  WDM  arrays . 

hybrid  and  monolithic  transmitters  using  a  100  km  of  Corning  SMF  28  single  mode  fibers 
and  a  detector  sensitivity  assessment  was  performed  by  comparing  our  results  at  a  BER  = 
1  x  10'9  with  and  without  the  fiber.  The  experimental  arrangement  is  shown  in  Fig.  4  and 
the  results  for  a  monolithic  Nortel  transmitter  are  shown  in  Fig.  5  a)  and  b.  There,  it  was 
found  that  the  penalty  could  be  negative,  of  order  1  dB,  for  some  transmitters,  indicating 
an  improvement  in  the  detector  sensitivity.  This  improvement  is  presumably  due  to  the 
sign  of  the  created  chirp,  which  might  create  some  pulse  compression,  therefore 
improving  the  signal -to-noise.  For  Ortel  transmitters,  a  penalty  of  order  1-2  dB  was 
typically  observed. 
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BER  Property  of  Nortel  Transmitter  98081 7-1-1  (No  Fiber  vs  100km  Fiber) 


Fig.  5  Transmission  experiment  evaluating  the  chirp  penalty  of  a  monolithic 
WDM  transmitter  array 

Crosstalk  Property  of  Norte!  Transmitter  98081 7-1  (100km  Fiber) 


Received  Power  (dBm) 

Figure  6.  Crosstalk  measurement  on  a  Nortel  monolithic  laser  array 


The  possibility  of  laser-to-laser  interactions,  cross-talk,  might  be  of  some 
concern  in  monolithically  integrated  devices.  Fig.  6  shows  the  results  of  a  transmission 
experiment  designed  to  evaluate  the  degree  of  cross  talk.  The  transmission  experiment 
measured  the  bit  error  rate  as  a  function  of  received  power,  through  100  km  of  fiber.  The 
experiment  was  first  carried  out  with  the  modulation  applied  only  to  one  of  the  lasers, 
laser  #1,  and  then  repeated  with  both  lasers  #1  and  #2  under  modulation.  The  presence  of 
cross  talk  would  normally  result  in  a  bit-error-rate-penalty,  i.e.  it  would  produce  a  shift 
between  the  two  lines.  The  experiment  illustrated  in  Fig.  6  shows  that  such  a  penalty  is 
indeed  very  small  in  the  monolithic  array  and  is  less  than  0.1  dB,  of  no  consequence  in 
most  transmission  experiments  of  interest.  This  degradation  mode  is  absent  in  hybrid- 
integrated  devices. 


Fig .  7  Long  term ,  1000  hrs ,  transmission  experiment  through  a  waveguide  router.  There 
were  no  errors  attributable  to  wavelength  or  power  instabilities. 


The  final  experiment.  Fig.  7,  illustrates  performance  of  our  lasers  in  a  test  WDM 
system.  In  this  experiment  the  laser  output  from  a  hybrid  array  was  first  combined  into  a 
single  fiber  and  connected  to  100  km  spool.  The  fiber  was  then  attached  to  a  5  x  5 
waveguide  grating  router  used  as  a  de-multiplexer.  The  router  is  capable  of  wavelength 
discrimination  greater  than  30  dB  and  a  offers  a  bandpass  of  ~0.5. nm.  Any  variation  in 
the  laser  output  outside  the  bandpass  of  the  router  would  result  in  catastrophic  error  rates. 
A  receiver  was  then  connected  to  the  center  output  waveguide  and  the  received  power 
was  adjusted  to  produce  a  bit-error-rate  of  IE-9.  The  laser  was  modulated  at  a  rate  of  2.5 
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Gbps.  The  setup  was  allowed  to  run  for  1000  hours  and  a  log  of  BER  measurements  was 
accumulated  under  computer  control.  The  highest  BER  recorded  in  that  period  was  6E-9. 

Summary 

We  have  investigated  wavelength  accuracy  and  stability  of  commercially  available 
WDM  transmitters.  Large  numbers  of  hybrid  and  monolithically  integrated  sources  were 
tested  as-delivered  from  manufacturers.  Transmitter  performance  was  tested  under  DC 
bias  and  modulation  at  data  rates  as  high  as  2.5Gbps.  The  wavelength  accuracy  of  hybrid 
devices  was  found  to  be  considerably  better  than  ±0.1  nm.  A  considerable  fraction  of 
devices  tested  showed  small  wavelength,  less  than  0.02  nm,  and  power,  less  than  0. 1  dB, 
excursions  due  to  digitization  noise  of  control  circuitry.  Transmission  experiments 
through  100  km  of  standard  single  mode  fiber  show  low  chirp  penalties  and,  in  case  of 
monolithic  arrays,  the  absence  of  cross-talk.  Long  term  transmission  experiments,  up  to 
1000  hrs,  through  wavelength-selective  elements  simulating  WDM  systems  do  not  show 
any  penalties  due  to  the  wavelength  or  power  instabilities.  The  commercial  devices  tested 
in  our  experiments  are  suitable  for  WDM  systems  research. 
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ABSTRACT 


A  method  of  propagation  functions  and  transfer  amplitudes  suitable  for  the  design  of  integrated 
optical  circuits  is  presented.  The  method  is  based  on  vectorial  formulation  of  electrodynamics:  the 
distributions  and  propagation  of  electromagnetic  fields  in  optical  circuits  is  described  by  equivalent 
surface  sources.  This  approach  permits  a  division  of  complex  optical  waveguide  structures  into  sets  of 
primitive  blocks  and  to  separately  calculate  the  transfer  function  and  the  transfer  amplitude  for  each 
block.  The  transfer  amplitude  of  the  entire  optical  system  is  represented  by  a  convolution  of  transfer 
amplitudes  of  its  primitive  blocks.  The  eigenvalues  and  eigenfunctions  of  arbitrary  waveguide  structure 
axe  obtained  in  the  WKB  approximation  and  compared  with  other  methods.  The  general  approach  is 
illustrated  with  the  transfer  amplitude  calculations  for  Dragone’s  star  coupler  and  router. 

Keywords:  integrated  optical  circuits,  waveguides,  star  couplers,  routers,  eigenmodes,  transfer 
amplitudes,  transfer  functions 


1  Introduction 


The  design  of  passive  waveguide  components  for  wavelength  division  multiplexing  and  integrated 
optical  circuits  requires  detailed  knowledge  of  light-guiding  parameters,  field  distributions,  transfer 
efficiencies  in  complex  device  structures.  Much  of  the  previous  work  done  in  this  field  concentrated  on 
evaluation  of  eigenvalues  and  eigenfunctions  for  specific  waveguide  structures.  The  finite-difference  and 
finite-element  methods,  the  effective  index  methods,  the  boundary  element  method,  the  WKB  method, 
numerous  numerical  methods,  etc.,  have  been  proposed.1"9  After  the  eigenvalues  and  eigenfunctions 
are  calculated,  their  propagation  through  an  optical  circuit  must  be  accounted  for.  Much  of  the  past 
effort  has  been  devoted  to  the  computation  of  field  distributions  and  propagation  characteristics  for 
devices  in  which  optical  energy  couples  between  the  input  and  output  ports.  The  best  known  examples 
are  the  work  of  Dragone,10-11  and  Henry  and  coworkers8-9  which  used  the  finite  Fourier  transform 
propagation  method.  Transform  techniques,  more  typical  of  electrical  engineering,  are  also  beginning 
to  be  applied  to  this  field.12 

The  previous  work  usually  considered  field  distributions  in  optical  waveguide  structures  in  terms 
of  scalar  electrodynamics.  In  this  work  we  develop  a  new  transfer-amplitude  method,  based  on  Love  s 


formulation13  of  vector  electrodynamics,  to  describe  field  distributions  and  transfer  functions  in  optical 
waveguide  systems.  This  approach  allows  us  to  divide  complex  optical  waveguide  structures  into  sets 
of  primitive  blocks  and  calculate  transfer  amplitudes  and  transfer  functions  for  each  block  separately. 
The  representation  of  a  complex  waveguide  structure  by  a  set  of  principal  surfaces  separating  primitive 
blocks  with  known  equivalent  surface  currents  and  propagation  functions  is  shown  to  be  a  powerful 
method.  It  provides  a  possibility  of  constructing  the  total  system  response  as  a  convolution  of  responses 
of  its  primitive  blocks.  The  exact  propagation  function,  including  diffraction  effects,  is  used  here  to 
describe  each  primitive  block,  instead  of  the  simple  ray  approximation.  We  apply  these  ideas  to  the 
evaluation  of  Dragone’s  star  coupler  and  router. 


2  Field  distributions  generated  by  equivalent  surface  sources 


Electromagnetic  field  solutions  of  Maxwell’s  equations  can  be  formulated  in  two  different,  but  equiv¬ 
alent  ways.  First,  the  time-harmonic  electromagnetic  fields  can  be  expressed  as  volume  integrals  over 
real  sources  which  are  vectors  of  real  current  densities.  On  the  other  hand,  the  solution  of  the  scalar 
Helmholtz  equation  can  be  represented  by  a  surface  integral  over  an  arbitrary  surface  outside  the  real 
source  region,  with  an  equivalent  surface  current  density.  The  latter  is  determined  by  the  field  gener¬ 
ated  by  the  real  source.  This  formulation,  called  Green’s  integral  theorem,  can  be  generalized  to  the 
case  of  vector  electrodynamics  and  the  result  is  known  as  Love’s  principle.1  In  this  Section,  we  present 
some  expressions  related  to  Love’s  principle  in  a  form  suitable  to  applications  in  optical  waveguides. 

Let  S'  (r)  =  0  be  the  equation  of  some  initial,  or  input,  surface  with  an  outer  normal  vector  »'. 
The  electric  field  of  the  TE  component  in  the  source-free  region,  i.e.,  outside  the  surface  S'  (r),  is 
determined  by  the  Helmholtz  equation  with  the  right  hand  side  expressing  the  boundary  condition, 


(V2  +  k%n2  (r))  E  (r)  =  -iu>n0Js  (r)  6  (S'  (r))  (1) 

where  the  equivalent  surface  current  density  is  given  by  Love’s  formula, 

Js  ( r )  =  2  n'  y.  H'  (r)  (2) 

Here  H'  (r)  is  the  magnetic  field  on  S'  produced  by  a  real  source;  u),  /r0,  ko  and  n  (r)  are  the  angular 
frequency,  the  permeability  of  vacuum,  the  wavenumber  of  fight  in  vacuum,  and  the  refractive  index, 
respectively;  r  is  a  radius  vector.  The  Dirac  delta-function  in  Eq.  (1)  allows  one  to  express  the  electric 
field  in  terms  of  Green’s  function,  G  (r,  r').  As  a  result,  the  volume  integration  reduces  to  integration 
over  the  surface  where  the  equivalent  surface  current  density,  Eq.  (2),  is  known, 


E  (r)  =  2iu>fi0  j  dS'G  (r,  r')  h'  x  H'  ( r ')  (3) 

S' 

A  similar  expression  can  be  generated  for  the  TM  component.  The  general  expression  for  the  field 
generated  by  an  equivalent  surface  source  can  be  then  simplified  in  cases  of  planar  waveguide  geometry 
when  a  standard  set  of  eigenfunctions  ip^z)  of  the  vertical  quantization  and  the  bare  propagation 
constant  are  given. 

Then,  the  Green’s  function  cam  be  reduced  to  a  two-dimensional  form, 


*'  =  0,1,2,  ... 


(4) 


where  p  and  p'  are  the  two-dimensional  radius  vectors  in  the  horizontal  plane,  and  v  is  the  mode 
number.  After  separating  the  longitudinal  and  vertical  coordinates,  the  electric  and  magnetic  fields  on 
the  input  surface  are  determined  by  the  eigenmodes  of  the  input  waveguide: 

E'  (r)  1 5'  =  e'E0(pS'  (Q  V>(z),  e'  =  s'  x  n'  (5) 


H'  (r)  |s.  =  —ip s>  (l)  +  Z'M(Z) 

(JfXQ 


(6) 


where  the  function  <Ps'{ 0  describes  the  field  distribution  on  the  principal  surface  S  (where  l  is  a 
transverse  coordinate);  the  mode  numbers  are  omitted  for  brevity;  z  is  a  unit  vector  along  the  vertical 
direction  and  e'  is  the  direction  of  the  electric  field  vector;  the  field  amplitude  Eo  is  an  arbitrary 
constant.  To  calculate  the  fields  outside  the  surface  S'  according  to  Eq.  (3),  one  needs  the  explicit 
form  of  the  equivalent  surface  current  density,  J*(r)  =  —e'%jj£Eo<Ps'(l)1l>(z)-  Thus,  the  general 
surface-integral  representation  (3)  can  be  reduced  to  a  linear  integral  over  the  transverse  coordinate  on 
the  surface  S': 


E  (r)  =  -2i(30Eorp  (z)  J  dl'e'G  (p,  p')  <pS'  (0  (7) 

S' 

where  V  is  the  arc  length  along  the  input  boundary  in  the  x-l  plane.  This  equation  describes  diffraction 
in  a  slab  region  between  the  input  and  output  principal  surfaces  in  a  given  primitive  block.  It  should 
be  emphasized  that  the  unit  vector  e!  appears  under  the  surface  integral  as  the  consequence  of  the 
vectorial  aspect  of  the  problem. 


3  Transfer  functions 


In  this  Section,  we  present  the  derivation  of  a  general  expression  for  the  amplitude  transfer  function 
(or,  in  brief,  the  transfer  amplitude)  based  on  Eq.  (7).  The  field  E(r)  is  generated  in  a  slab  region 
by  an  effective  surface  source  on  S'.  It  propagates  from  an  input  surface  S'  to  an  output  surface 
S  according  to  Huygens’  principle  and  enters  the  system  of  output  waveguides.  Let’s  suppose  the 
total  eigenfunction  of  the  output  waveguide  system  contains  a  product  of  a  new  unit  vector  e  and  an 
eigenfunction  describing  vertical  quantization  ip(z).  This  representation  separates  the  vector  properties 
and  the  vertical  coordinate  dependence.  The  remaining  spatial  dependence,  related  to  the  transverse 
direction,  is  partly  described  by  an  eigenfunction  <Ps(l)-  Besides,  the  electric  field  at  the  output  surface 
as  given  by  a  product  of  the  vertical  and  transverse  eigenfunctions  is  reduced  by  one  additional  factor, 
i.e.,  the  transfer  amplitude,  describing  diffraction  in  the  primitive  block. 

A  natural  way  to  define  a  dimensionless  transfer  amplitude  for  a  primitive  block  consisting  of  surfaces 
S  and  S'  and  a  slab  region,  is  to  consider  the  functional  and  vectorial  projections  of  the  electric  field 
E(r),  as  defined  by  Eq.  (7),  onto  the  product  of  eigenfunctions  of  the  output  waveguides,  divided  by 
a  field  amplitude  Eq: 


t{s,s')  =  Y0JdStps( woo  s  •£(»•) 


(8) 


(9) 


=  -2»|80  J  dl  J  dl'e-eYsWGfop'M') 

S  S' 

To  get  the  second  line  in  this  equation,  we  used  the  orthonormality  properties  of  the  eigenfunctions  ip 
and  <p.  With  the  definition  (8),  the  electric  and  magnetic  fields  related  to  the  output  surface  5  can  be 
represented  as  follows, 


E  (r)|s  =  eEot  (S,  S’)  <ps  (1 l )  ip  (*)  (10) 


H  (r)|5  =  ^-t  (S',  S')  <ps  (J)  +  *M>  (*)  (n) 

These  expressions  for  outgoing  fields  on  the  surface  S  are  very  similar  to  expressions  for  the  incoming 
fields  on  the  surface  S',  Eqs.  (5)  and  (6),  except  for  the  presence  of  the  amplitude  transfer  factor. 

The  normalization  condition  which  is  necessary  for  the  construction  of  a  transfer  function  requires 
the  knowledge  of  the  input  Poynting  vector  on  the  surface  S',  The  Poynting  vector  of  the 

outgoing  field  can  be  easily  obtained  on  the  basis  of  Eqs.  (10)  and  (11).  After  integration  over  the 
surface  S,  we  get 


Ptotal  =  h^-  \Eo\2  |t  (5,  S")|2  =  Ptotai  1 1  {S,  S') I2  (12) 

Thus,  the  quantity  1 1  (5,  S")|2,  where  t  is  defined  by  Eq.  (8),  represents  a  power  (or  intensity)  transfer 
function  for  a  given  primitive  block, 


T(S,S')  =  \t{S,S')\2  (13) 

The  total  transfer  amplitude  for  the  entire  optical  system  is  represented  by  a  convolution  of  multiple 
transfer  amplitudes  of  primitive  blocks,  with  summation  over  all  intermediate  principal  surfaces  (Fig. 

1): 


t{Sout,Sin)  =  52  t(SouUSn)...t(S2,S1)t(S1,Sin)  (14) 

Si  ,S2,-.mStx 

The  corresponding  transfer  function  for  the  entire  optical  system  is  defined  as. 


T  (Sent,  Sin)  =  1 1  (S^t,  Sin) |2  (15) 

At  this  point,  it  is  convenient  to  introduce  a  normalized  transverse  eigenfunction  of  the  diffracted 
field,  according  to  Eq.  (9): 

< p(l )  =  —  2i/?0  J  dlf  e  •  e'G(p , 

S' 

The  concept  of  a  normalized  transverse  eigenfunction  of  the  diffracted  field  is  very  useful  because  (i) 
the  function  <p{l)  is  defined  uniquely  by  Eq.  (16)  everywhere,  (ii)  it  can  be  easily  generalized  to  the  case 


of  several  primitive  blocks,  and  (iii)  the  transfer  amplitude,  Eq.  (9),  can  be  represented  as  a  complex 
overlap  integral  of  two  normalized  eigenfunctions  ip(l)  and  <p$(l)  on  the  principal  surface  5,  i.e., 

t(S,S*)  =j  dl<fi*s(l)  ¥>(*)  (17) 

s 

These  two  eigenfunctions  describe  the  incoming  and  outgoing  fields  on  the  principal  surface  S. 

In  calculating  the  properties  of  integrated  optical  devices,  it  is  easy  to  see  that  the  eigenfunctions 
ips(l)  and  ip{l)  should  closely  match  each  other  on  every  principal  or  output  surface  S  in  order  to 
maximize  efficiency  or  the  transfer  function,  Eq.  (13)  or  (15),  of  a  particular  design.  In  other  words, 
each  pair  of  these  two  functions  (for  a  given  S)  has  to  be  chosen  commensurable. 


4  WKB  approximation  for  optical  waveguides 


The  spatial  dependence  in  the  transverse  direction,  described  above  by  an  eigenfunction  y>s(i)»  is 
generally  not  known.  This  function  is  a  solution  of  a  second  order  differential  equation  and  it  can 
be,  in  principle,  calculated  exactly.  It  can  be  also  modeled  for  specific  waveguide  designs  by  assuming 
a  cosine  or  more  complex  periodic  refractive  index  or  aperture  profiles.10  This  Section  describes  a 
method  of  calculating  the  transverse  eigenfunction  <ps(l)  which  is  based  on  WKB  approximation  (see, 
for  example,2’14).  It  is  assumed  that  the  vertical  quantization  is  performed  first  and  that  the  index 
profile  in  the  transverse  direction  is  a  smooth  function.  We  are  thus  able  to  use  the  known  vertical 
eigenfunctions  and  obtain  <ps(l)  from  the  standard  WKB  approximation.  The  electromagnetic  field 
dependence  on  the  longitudinal  coordinate  (x)  is  determined  by  an  effective  propagation  constant  /3, 
to  be  determined  below. 

The  WKB  approach  developed  here  can  be  formulated  in  the  following  manner.  For  a  waveguide 
system  of  arbitrary  geometry  in  the  x-l  plane  (the  local  Cartesian  coordinates  axe  x,  l,  z)  and  a 
rectangular  index  step  in  the  z  direction,  a  general  eigenmode  can  be  described  as 


E(r)=eE0eil3x<p{l)ip(z)  (18) 

where  /?  and  /30  (/30  >  j9)  are  the  effective  and  bare  propagation  constants,  respectively,  e  is  a  unit 
vector,  and  E0  is  a  constant  initial  amplitude.  We  look  for  a  transverse  eigenfunction  as  given  by  the 
usual  WKB  expression,  tp  (Z)  =  Re  (e*fc°s^).  This  eigenfunction  describes  the  waveguide  region  between 
the  principal  surfaces,  so  that  the  surface  index  (5  or  S')  is  inappropriate  here.  After  separation  of  the 
longitudinal  coordinate,  the  product  of  ip  and  ip  in  Eq.  (18)  should  satisfy 

{w + £ + (*°v  z)  ~ p2)) (v?  (0  * {z)) = 0  (19) 

In  addition  ip  must  satisfy  the  condition  of  vertical  quantization: 

( +  kon*  (*)  -  *)  ^  (z)  =  0  (20) 


We  assume  a  general  form  of  the  refractive  index  profile  in  the  transverse  direction  as  given  by 


where  «*(/)  represents  a  single-well  (or  single-waveguide)  or  multiple-well  structure  in  the  transverse 
direction  l.  After  substituting  the  general  form  of  a  solution  into  the  system  of  equations  (19)  and  (20), 
the  basic  WKB  equation  for  the  eigenfunction  <p{l)  is  obtained: 


(|U^p2(o)  ^)=0  (22) 

where  p{l)  is  the  analogue  of  the  classical  momentum  of  a  particle  in  the  quantum  mechanical  WKB 
approximation, 


p(‘)  =  p^-M -<('))]*  (23) 

The  Hiffiran™.  between  the  formulation  of  the  WKB  approximation  presented  here  for  optical 
waveguides  and  the  standard  formulation,  as  described  for  example  in,2  is  that  the  condition  of  vertical 
quantization  is  included  explicitly,  Eq.  (20).  This  results  in  the  appearance  of  a  difference  between 
propagation  constants  in  p(Z).  After  series  expansion  in  the  characteristic  parameter  ,  the  WKB 
eigenfunction  can  be  expressed  as 


/  ±iko  f  dxP(i)  \ 

¥’(0  =  7&)Re(e  11  )’  h<l<h  (24) 

where  h  and  l2  are  the  adjacent  classical  turning  points  determined  by  the  roots  of  the  function  p(Z); 
Q  is  a  normalization  constant.  In  classically  forbidden  regions,  the  momentum  becomes  imaginary  and 
the  oscillating  exponential  in  (24)  reduces  to  the  decaying  exponential.  As  usual,  the  WKB  solution 
(24)  is  valid  only  far  away  from  the  turning  points.  If  this  solution  is  formally  applied  to  the  regions 
near  the  turning  points,  artificial  singularities  appear  due  to  the  fact  that  the  classical  momentum  p(Z) 
equals  zero  at  each  turning  point.  These  singularities  are  integrable  and  do  not  affect  the  eigenfunction. 

Lets  consider  a  case  in  which  the  function  ntr(Z)  in  Eq.  (21)  describes  the  index  profile  of  a  symmetric 
two-waveguide  structure  with  two  classically  allowed  and  three  forbidden  regions.  The  turning  points, 
lu  are  the  boundaries  between  these  regions.  The  comparison  of  solutions  in  the  allowed  and  forbidden 
regions  gives  rise  to  quantization  conditions  for  symmetric  and  asymmetric  eigenfunctions,  respectively. 


*2 


k0  f  p(l )  dl  =  pit  +  arctan(2e“)  (symmetric  case) 


(25) 


(2 

fc„  f  p(l)  dl  =  (p  +  l)7r  -  arctan(2e“)  (asymmetric  case) 


(26) 


where  p  is  zero  or  a  positive  integer  and  the  parameter  a  =  k0  J  |p(Z)|  dl.  Eqs.  (25)  and  (26)  represent 

a  generalization  of  the  Bohr-Sommerfeld  quantization  condition  to  the  case  of  a  double-well  structure. 
It  is  instructive  to  see  how  this  more  general  solution  reduces  to  the  proper  Bohr-Sommerfeld  condition 
in  the  case  when  the  separation  between  two  quantum  wells  tends  to  infinity.  In  this  case,  a  tends 


to  infinity  (the  limit  of  an  infinitely  long  potential  barrier),  arctan(2ea)  approximately  equals  §,  and 
both  equations  (25)  and  (26)  reduce  to  a  single  equation  of  the  following  form 


*  /  \ 

ko  J p(l)  dl=  \  p+  -Jn  (27) 

h 

with  double  degeneracy  for  each  eigenvalue  for  a  given  If  we  neglect  this  degeneracy,  Eq.  (27) 
becomes  the  standard  Bohr  Sommerfeld  condition.  Consequently,  Eqs.  (25),  (26)),  and  (27)  provide 
the  basis  for  numerical  calculations  of  the  effective  propagation  constant  0  for  optical  single-well  and 
double-well  waveguide  structures. 

The  effective  index  method  is  based  on  successive  application  of  quantization  in  the  vertical  and 
transverse  directions.  This  method  can  be  applied  to  the  refractive  index  profile  with  abrupt  boundaries 
in  both  of  these  directions.  The  method  of  tins  paper  is  more  general  because  it  is  applicable,  with  good 
accuracy,  to  any  index  profile  in  the  transverse  direction  including  the  case  of  a  smooth  profile.  The 
procedure  of  vertical  quantization  within  the  present  approach  is  exactly  the  same  as  in  the  effective 
index  method.  The  effective  propagation  constant  of  a  single-waveguide  structure  is  presented  in  Fig.  2 
as  a  function  of  the  half-width  of  the  waveguide,  a,  and  compared  with  that  obtained  using  the  effective 
index  method.  The  results  of  numerical  calculations  of  the  effective  refractive  index  for  the  cases  of 
hyperbolic  cosine  profile  and  the  almost  rectangular  waveguide  profile  are  also  presented  in  Fig.  2. 
Figure  2  shows  very  good  agreement  between  the  WKB  approximation  for  the  hyperbolic  cosine  index 
profile  (curve  1)  and  the  result  of  the  effective  index  method  for  a  rectangular  profile  (curve  3).  The 
difference  between  the  two  calculations  is  only  about  10~3.  To  extend  this  comparison  further  we  also 
show  the  exact  numerical  calculation  (curve  2)  for  the  profile  used  in  the  WKB  approximation.  It  is 
easy  to  see  that  curves  1  and  2  are  very  close  for  a  >  1  /ttm.  In  addition,  we  present  results  of  two 
calculations  for  an  almost  rectangular  profile  done  by  the  WKB  method  and  the  numerical  calculation 
(curves  4  and  5,  respectively),  to  be  compared  with  the  effective  index  method  results  (curve  3).  The 
effective  index  method  is  in  good  agreement  with  the  numerical  calculations  (curve  5),  as  expected,  and 
the  difference  between  curves  3  and  5  is  almost  negligible  for  a  >  2  fixn.  The  WKB  approximation  for 
the  rectangular  profile  (curve  4)  also  works  well  for  a  >  2  /im,  where  it  differs  by  about  0.01%  from  the 
exact  numerical  calculation.  However,  for  this  abrupt  index  profile  the  WKB  approximation  cannot  be 
used  for  a  <  2  (im.  The  advantage  of  the  WKB  method  is  the  possibility  of  calculating  eigenfunctions 
corresponding  not  only  to  single  waveguides  with  the  effective  propagation  constant  0  (which  can  be 
obtained  from  Fig.  2  with  0  -  fco«e//).  but  of  extending  such  calculations  to  more  general  transverse 
waveguide  structures.  This  is  in  contrast  with  the  effective  index  method. 

In  Figure  3  we  illustrate  a  procedure  for  eliminating  singularities  at  the  classical  turning  points  of 
a  single-waveguide  structure  and  compare  the  WKB-eigenfunctions  with  the  eigenfunctions  calculated 
numerically.  Figure  4  shows  the  modified  WKB  eigenfunction  (curve  1)  and  the  corresponding  exact 
eigenfunction  (curve  2)  for  the  case  of  a  double-well  waveguide  structure,  where  2d  is  the  distance 
between  individual  waveguide  centers  {d  =  3  fJ.m) . 


5  Applications  to  star  couplers  and  routers 


In  this  Section  we  consider  propagation  of  eigenmodes  and  transverse  eigenfunctions  using  the 
transfer  amplitude  defined  in  Eq.  (17).  It  is  well-known  that  diffraction  of  guided  waves  in  a  star 
coupler  with  circular  (cylindrical)  principal  surfaces  and  the  ratio  d/R  1  (where  d  is  the  distance 
between  waveguides  in  a  periodic  structure  in  the  transverse  direction,  and  R  is  the  characteristic 
radius  of  principal  surfaces  bounding  the  slab  region)  can  be  described  with  high  accuracy  within  the 


Fourier-optics  approximation.10*15  The  calculations  of  the  transverse  eigenfunction  <p{l)  of  the  diffracted 
field  and  the  transfer  amplitude  t  (S,  S')  for  Dragone’s  star  coupler  and  router,  i.e.,  two  star  couplers 
connected  through  a  multi-waveguide  grating  structure,  are  considered  below  in  this  approximation. 

In  the  case  of  a  star  coupler,  the  basic  Eqs.  (16)  and  (17)  can  be  reformulated  as. 


t  (5,  S')~J  d®  d>s(a)  <Ka)  (28) 


<H«)  = 


da'e-^0R  aa'<t>S‘(<x') 


(29) 


where  the  integral  in  Eq.  (29)  is  a  typical  Fourier-optics  integral.  Here  we  use  angular  variables  a 
and  o'  (a,  o'  <  1,  e  •  e'  «  1)  instead  of  a  pair  of  two-dimensional  radius-vectors  p  and  p ,  and  the 
eigenfunctions  are  normalized  by  the  angle  a  instead  of  the  arc  length.  If  the  input  into  the  slab  region 
is  introduced  by  a  single  waveguide  with  the  ordinal  number  m  on  S',  the  transverse  eigenfunction 
(29)  of  the  diffracted  field  on  S  can  be  represented  in  general  by  a  Bloch  function  of  the  order  m: 
<b  (a)  =  e~im0oda  <£„(<*).  The  angular  period  of  these  Bloch  oscillations  equals  27r/(|m|/?0d)  and  is 
much  smaller,  for  |m|  >  1,  than  the  characteristic  spread  of  the  diffracted  field  described  by  the  function 
^0(a)  (this  spread  is  of  the  order  4x /(/?0d)).  On  the  other  hand,  this  period  is  much  greater  than  the 
natural  angular  period  of  the  waveguiding  system,  d/R. 


When  the  eigenfunctions  <f>0{a)  and  <f>s(a)  Me  hnown»  the  integral  in  Eq.  (28)  can  be  easily  eval¬ 
uated.  The  calculated  transfer  functions  T  =  \tf  are  presented  in  Fig.  5  for  different  values  of  N 
Transfer  functions  are  sensitive  to  the  size,  the  total  number,  and  the  periodicity  of  the  input  an 
output  waveguides.  This  allows  for  easy  optimization  of  the  device.  The  calculated  values  for  the  star 
coupler’s  efficiency  are  realistic  and  relatively  high,  reaching  a  maximum  of  ~  80%.  That  is  much 
higher  than  the  values  estimated  previously.10 

In  the  case  of  Dragone’s  router,  from  the  purely  geometrical  consideration,  one  can  establish  that 
the  transmission  through  a  multi-waveguide  system,  from  Si  to  S[  can  be  described  by  a  special  type 
of  the  transfer  amplitude,  called  here  W^ai),  which  is  related  to  the  grating  geometry  and  the 
detailed  shape  of  the  waveguide  openings  at  the  two  surfaces.  When  both  star  couplers  are  of  the  same 
shape  and  size,  the  transfer  amplitude  W  is  proportional  to  the  delta  function  of  its  arguments  and  we 
get  the  following  symmetric  form: 


W(aj,ai)  =  v)(ai)  £(ai  —  <*i)  (30) 

where  the  single-parameter  phase  factor  w{ax)  describes  phase  shifts  between  adjacent  waveguides 
arising  from  their  relative  displacements  and  different  optical  paths  between  S\  and  Si .  Using  Eq.  (30) 
the  general  formulas  describing  double  diffraction  in  two  star  couplers  can  be  represented  as. 


tiSouuSin)^  J  da  <psout(a)  4>(a)  (31) 


da1  K  {a,  a1)  <pSin  (a') 


(32) 


(33) 


K(a,a')  =  (—2i/30R)2  J  da i  G(a,a!)  w(a x)  G (aq,a') 

Si 

The  remaining  double  integral  in  Eqs.  (32)  and  (33)  is  quite  simple  in  the  Fourier-optics  approximation. 
The  grating  produces  a  periodic  signal  on  the  output  surface  Sout •  In  this  case,  the  characteristic 
function  w  can  be  represented  by  a  sum  of  discrete  contributions  from  each  waveguide  in  the  grating. 


w(c*i)  =  52  eis/3A  6(ai  -  sad)  (34) 

S 

where  A  is  the  length  difference  between  adjacent  waveguides,  the  product  of  /?A  is  the  optical  path 
difference  between  waveguides,  and  s  =  0,±1,  ±2, ....  Integrating  over  the  angular  variables  and  using 
the  well-known  formula  related  to  Floquet’s  theorem  we  get: 


K  (a ,  a1 )  —  —ie2iP°R  52  $(<*  +  «'- 

p=— OO 


I 3  A  2tt p 

Pod  Pod’ 


The  grating  thus  generates  a  periodic  delta  function  in  the  transfer  amplitude. 


(35) 


Eqs.  (32)  and  (35)  have  important  applications  to  wavelength  multiplexing  and  demultiplexing  in 
integrated  optical  circuits.  Suppose  the  input  to  the  router  consists  of  four  closely  spaced  and  equidis¬ 
tant  wavelengths  in  a  single-waveguide,  with  four  corresponding  values  of  the  propagation  constant 
p.  This  set  of  wavelengths  can  be  represented  by  A„  =  Ao(l  +  ne)  where  n  =  0, 1,2,3,  and  Ao  and 
pQ  the  initial  values  for  the  (demultiplexing  problem.  To  resolve  any  two  adjacent  wavelengths, 
the  condition  p0  =  on  the  last  term  in  the  argument  of  the  delta  function  in  Eq.  (35)  must  be 

satisfied,  where  p0  is  the  order  of  the  diffraction  maximum.  Then  the  a'  in  the  delta  function  in  Eq. 
(35)  takes  the  following  values: 


a'  =  -a-nad  +  (0-^)  (36) 

where  n  =  0, 1, 2, 3.  The  light  of  any  wavelength  can  be  thus  send  to  any  output  waveguide,  depending 
only  on  the  input  waveguide  used  and  the  parameter  A.  The  phase  shift  s*  in  the  argument  of  (36) 
can  be  partly  or  entirely  compensated  for  by  choosing  the  value  of  the  characteristic  path  difference 
A.  For  example,  complete  compensation  is  achieved  for  A  =  ^  ^f-d.  In  this  situation,  the  light  of 
wavelength  Ao  entering  the  central  input  waveguide  is  directed  to  the  central  output  waveguide.  The 
same  light  entering  the  waveguide  m  is  directed  to  the  output  waveguide  -m.  Light  of  the  wavelength 
Ai  is  displaced  by  one  angular  period  of  the  grating,  etc. 

The  different  wavelength  contributions  are  totally  separated  (demultiplexed)  in  space  (angular) 
variable  and 


3 

<j>{a)  ~  52  ^sin  (~a  -  nad )  (37) 

n=0 

where  c„  is  the  amplitude  coefficient  of  a  contribution  with  An  in  the  input  signal.  We  conclude  that 
the  input  signal  is  totally  demultiplexed  on  the  output  surface  of  Dragone’s  router,  i.e.  each  wavelength 
enters  its  own  output  waveguide.  When  the  phase  shift  in  Eq.  (36)  is  not  completely  compensated  for, 
an  additional  phase  shift  should  be  entered  into  Eq.  (37).  This  might  be  useful  when  it  is  desirable 
to  send  any  wavelength  to  any  particular  output  waveguide.  The  typical  values  A0  =  1.5500  pm  and 
AA  =  0.0032  pm  correspond  to  the  diffraction  maximum  of  the  order  po  =  106  in  Eq.  (37). 


An  example  of  wavelength  demultiplexing  is  shown  in  Fig.  6  which  plots  the  diffracted  intensities 
for  two  wavelengths  as  a  function  of  the  output  angle.  The  resolving  capability  of  the  device  depends 
on  the  total  number  of  waveguides  in  the  grating  section.  Figure  6  shows  the  response  of  the  router 
with  19  waveguides  in  the  grating. 


6  Summary 


In  conclusion,  the  method  of  modeling  complex  optical  systems  developed  in  this  work  provides  a 
qualitative  explanation  for  the  propagation  of  eigenmodes  in  integrated  optics  devices  and  also  a  very 
powerful  and  flexible  method  for  numerical  calculations  of  diffracted  field  distributions,  efficiencies,  etc., 
needed  for  complete  characterization  of  such  devices.  The  proposed  design  technique  consists  of  the 
following  steps:  (a)  dividing  a  complex  optical  waveguide  system  into  a  set  of  primitive  blocks  separated 
by  principal  optical  surfaces,  and  the  determination  of  boundary  values  of  electromagnetic  fields  and 
their  polarization,  (b)  expressing  the  diffracted  fields  through  surface  integrals,  generated  according  to 
Love’s  field  equivalence  principle,  by  the  equivalent  surface  current  densities,  (c)  computing  the  transfer 
amplitude  and  transfer  function  for  each  block  separately,  (d)  calculating  the  transfer  amplitude  of  the 
entire  optical  system  as  a  convolution  integral  of  the  transfer  amplitudes  of  its  primitive  blocks,  (e) 
evaluating  the  efficiency  of  an  optical  device,  i.e.,  its  total  transfer  function,  as  a  function  of  the  size 
and  shape  of  apertures  on  principal  surfaces  or  the  number  of  waveguides  in  the  grating. 
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Figure  Captions 

Figure  1:  Schematic  representation  of  an  integrated  optical  structure  with  its  input,  output  and 
intermediate  principal  surfaces  separating  the  primitive  blocks,  representative,  for  instance,  of  a  Drag- 
one’s  router.  The  Nij  stands  for  a  number  of  waveguides  between  the  ith  and  jth  principal  surfaces. 

Figure  2:  The  effective  refractive  index  (ne//  =  P/ko)  of  a  single-well  waveguide  plotted  as  a 
function  of  its  half-width  a,  as  calculated  using  the  WKB  approximation,  the  effective  index  method, 
and  the  exact  numerical  calculations:  1  -  the  case  of  hyperbolic  cosine  profile  calculated  by  the  WKB 
approximation  (solid  line)  and  2  -  the  numerical  calculation  (A);  3  -  the  case  of  a  rectangular  profile 
with  abrupt  boundaries  calculated  by  the  effective  index  method  (dashed  fine);  4  -  the  case  of  an  almost 
rectangular  profile  calculated  by  the  WKB  approximation  (□)  and  5  -  the  numerical  calculations  (■). 
The  horizontal  asymptotic  line  represents  the  value  of  the  effective  index  in  the  slab  region  (/3q/^o)- 
All  the  calculations  use  parameters  characteristic  of  semiconductor  waveguides. 

Figure  3:  The  TE  eigenmodes  of  a  single  waveguide  with  the  hyperbolic  cosine  profile  calculated 
by  (1)  direct  application  of  the  WKB  approximation,  (2)  after  elimination  of  singularities  at  classical 
turning  points,  and  (3)  by  numerical  calculations  . 

Figure  4:  The  TE  eigenmodes  for  coupled  waveguides,  (1)  after  elimination  of  singularities  at  the 
four  classical  truning  points  ,  and  (2)  by  numerical  calculations  . 

Figure  5:  The  efficiency,  T  (5,  Sf ) ,  of  Dragone’s  star  coupler  plotted  as  a  function  of  the  number  iV 
of  output  waveguides  in  the  case  of  a  single- waveguide  input:  1  and  2  -  the  case  of  a  merging  waveguide 
system  on  S  with  a  =  d/2  and  a  =  d/4,  respectively;  3  -  the  case  of  a  periodic  cos2(^)  aperture 
profile  on  S. 

Figure  6:  The  resolving  capability  of  two  Dragone’s  couplers  with  19  waveguides  in  the  grating,  in 
the  case  of  a  single-waveguide  input:  the  diffracted  intensities  for  two  wavelengths  as  a  function  of  the 
output  angle  according  to  Eqs.  (35)  and  (37). 
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Propagation  of  Eigenmodes  and  Transfer 
Amplitudes  in  Optical  Waveguide  Structures 

V.  A.  Mashkov  and  H.  Temkin,  Fellow \  IEEE 


Abstract — A  method  of  transfer  amplitudes  suitable  for  mod¬ 
eling  and  simulation  in  integrated  optical  circuits  is  presented. 
The  method  is  based  on  vectorial  formulation  of  electrodynamics: 
distributions  and  propagation  of  electromagnetic  fields  in  optical 
circuits  are  described  by  equivalent  surface  sources.  This  ap¬ 
proach  permits  a  division  of  complex  optical  waveguide  structures 
into  sets  of  primitive  blocks  and  to  separately  calculate  the 
transfer  amplitude  for  each  block.  The  transfer  amplitude  of  the 
entire  optical  system  is  represented  by  a  convolution  of  transfer 
amplitudes  of  its  primitive  blocks.  Diffraction  in  slab  regions 
is  taken  into  account  using  known  propagation  functions.  It  is 
shown  that  the  concept  of  a  normalized  transverse  eigenfunction 
of  the  diffracted  field  is  very  useful  in  calculations  of  transfer 
amplitudes.  The  crucial  role  of  transverse  Bloch  modes  in  prop¬ 
agation  through  multiwaveguide  structures  is  emphasized.  With 
this  method,  the  eigenvalues  and  eigenfunctions  of  an  arbitrary 
waveguide  structure  can  be  obtained  with  high  accuracy.  The 
general  approach  is  illustrated  with  the  transfer  amplitude  and 
efficiency  calculations  for  a  star  coupler  and  a  waveguide  grating 
router. 

Index  Terms — Eigenmodes,  integrated  optical  circuits,  star  cou¬ 
plers,  transfer  amplitudes,  waveguide  grating  routers,  waveg¬ 
uides,  wavelength  division  multiplexing  and  demultiplexing. 


I.  Introduction 

THE  DESIGN  OF  passive  waveguide  components  for 
wavelength  division  multiplexing  and  integrated  optical 
circuits  requires  detailed  knowledge  of  light-guiding  parame¬ 
ters,  field  distributions,  transfer  efficiencies,  etc.,  in  complex 
device  structures.  Much  of  the  previous  work  done  in  this 
held  concentrated  on  two  aspects  of  this  problem.  First, 
several  powerful  methods  have  been  developed  to  evaluate 
eigenvalues  and  eigenfunctions  for  specific  waveguide  struc¬ 
tures.  These  are  the  finite-difference  or  multigrain  method, 
the  finite-element  or  variational  method,  the  effective  index 
methods,  the  boundary  element  or  integral-equation  method, 
and  a  specific  matrix  implementation  of  it  known  as  the 
Galerkin  method,  the  WKB  method,  numerous  numerical 
methods,  etc.  (see  [1|-[I3]  and  references  therein).  Second, 
after  calculating  the  eigenvalues  and  eigenfunctions,  their 
propagation  through  an  optical  circuit  must  be  accounted 
for.  The  propagation  aspect  makes  the  calculation  of  field 
distributions  a  complicated  and  time-consuming  task.  Much 
of  the  past  effort  has  been  devoted  to  the  computation  of  field 
distributions  and  propagation  characteristics  for  optoelectronic 
devices  in  which  optical  energy  couples  between  the  input 
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and  output  ports.  The  best  known  examples  are  the  work  of 
Dragone  [14],  [15]  and  Henry  et  al  [11],  [12],  which  used 
Fourier  transform  propagation  methods.  In  multiwaveguide 
structures  with  gradual  transitions,  such  as  waveguide  merging 
or  tapering,  propagation  characteristics  of  optical  fields  can 
be  described  using  Bloch  modes  [16].  Transform  techniques, 
more  typical  of  microwave  engineering,  are  also  beginning  to 
be  applied  to  this  field  [17]. 

The  previous  work  usually  considered  field  distributions 
in  optical  waveguide  structures-in  terms  of  scalar  electrody¬ 
namics.  However,  the  pioneering  work  on  radar,  antennas, 
and  microwave  waveguides  has  demonstrated  the  general 
usefulness  of  describing  electromagnetic  fields  in  terms  of 
vector  electrodynamics  [2].  The  significance  of  this  comes 
from  the  fact  that  optical  waveguide  systems  are  described 
not  only  by  the  wave  equation  for  the  fields,  but  also  by  the 
boundary  conditions  which  are  vectorial  in  nature. 

In  this  paper,  we  develop  a  transfer-amplitude  method  to 
describe  field  distributions  and  transfer  functions  in  optical 
waveguide  systems  based  on  vector  electrodynamics.  The 
starting  point  of  the  analysis  is  Huygens’  picture  of  prop¬ 
agation  of  time-harmonic  electromagnetic  fields  in  complex 
waveguide  structures  as  expressed  by  Love’s  formulation  1 1 8 1 
of  the  vectorial  Green’s  theorem.  This  approach  allows  us  to 
divide  a  complex  optical  waveguide  structure  into  a  set  of 
primitive  blocks  and  to  calculate  the  transfer  amplitude  for 
each  block  separately.  We  demonstrate  that  Love's  expressions 
for  equivalent  surface  current  densities  lead  to  correct  results 
for  the  field  distributions  in  primitive  blocks.  It  is  then 
possible  to  calculate  partial  and  total  transfer  amplitudes  and 
corresponding  device  efficiencies  with  high  accuracy. 

Section  II  describes  field  distributions  generated  by  equiva¬ 
lent  surface  sources  based  on  Love’s  principle.  We  show  how 
to  extract  the  dependence  of  the  electromagnetic  field  on  the 
vertical  coordinate,  the  confinement  direction,  and  to  obtain 
two-dimensional  (2-D)  field  distributions. 

General  expressions  for  transfer  amplitudes  of  primitive 
blocks  and  of  the  entire  optical  system  are  presented  in 
Section  III.  A  waveguide  structure  is  represented  by  a  set 
of  principal  surfaces  separating  primitive  blocks,  each  with 
known  equivalent  surface  currents  and  propagation  functions. 
The  total  system  response  is  then  constructed  as  a  convolution 
of  responses  of  its  primitive  blocks. 

The  effects  of  lateral  confinement  and  transverse  eigenfunc¬ 
tions  for  single  and  multiple  waveguides  are  considered  in 
Section  IV,  where  simulations  of  single,  double,  and  periodic 
waveguide  structures  are  presented. 
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Sections  V  and  VI  present  the  examples  of  applications  in 
which  multiple  diffraction  effects  are  important.  The  general 
ideas  are  applied  to  the  analysis  and  efficiency  evaluation  of  a 
star  coupler  and  Dragone’s  router  (two  star  couplers  connected 
by  a  waveguide  grating)  with  emphasis  bn  the  Bloch-mode 
nature  of  diffracted  fields. 

II.  Field  Distributions  Generated 
by  Equivalent  Surface  Sources 

Electromagnetic  field  solutions  of  Maxwell’s  equations  can 
be  formulated  in  two  different  but  equivalent  ways.  First, 
time-harmonic  electromagnetic  fields  Can  be  expressed  as 
volume  integrals  over  real  sources  which  are  vectors  of  real 
conservative  current  densities.  Second,  the  solution  can  be 
represented  by  a  surface  integral  over  an  arbitrary  surface 
outside  the  real  source  region,  with  an  equivalent  surface 
current  density  [19].  The  latter  is  determined  by  the  fields 
generated  by  the  real  source.  This  formulation  is  known  as 
Love’s  principle  of  vector  electrodynamics  [2],  [19],  and  it 
can  be  reduced  to  Green’s  integral  theorem  in  case  of  scalar 
electrodynamics. 

Love’s  field-equivalence  principle  uses  surface  currents  to 
describe  optical  fields  in  the  source-free  region  [21.  The  optical 
fields  generated  by  any  real  source  and  the  optical  fields 
generated  by  the  corresponding  equivalent  surface  sources 
are  totally  equivalent.  Due  to  the  boundary  conditions,  the 
presence  of  the  surface  itself  decouples  the  total  Optical  field 
into  vector  components,  subject  to  Fresnel’s  electric  and 
magnetic  refraction  at  the  surface  [19].  This  decomposition 
is  unique  and  it  is  not  subject  to  any  approximations.  In 
case  of  a  curved  or  uneven  surface,  this  decomposition  is 
done  at  each  point  on  the  surface,  i.e.,  it  is  local.  The  two 
resulting  polarization  components  are  thus  described,  outside 
the  surface,  by  two  independent  surface  currents  defined  at 
each  point  on  the  surface.  Love’s  principle  provides  rules  for 
the  evaluation  of  equivalent  surface  sources  when  the  optical 
fields  on  a  given  principal  surface  are  known  [2],  118]. 

We  formulate  Love’s  principle  by  analogy  with  the  treat¬ 
ment  of  bulk  current  densities  in  Maxwell’s  equations.  Let 
S'( r)  =  0  be  the  equation  of  some  initial,  or  input,  surface 
with  an  outer  unit  normal  vector  n'  (Fig.  1).  Due  to  the  surface 
decomposition  described  above,  there  are  two  surface  currents 
and  two  distinct  modes.  One  of  the  modes  is  very  similar  to  a 
TE  mode  in  a  slab  waveguide,  while  the  other  one  is  similar  to 
a  TM  mode.  We  use  this  similarity  to  simplify  the  terminology 
throughout  the  paper:  we  refer  to  TE-like  modes  simply  as  TE 
modes,  and  TM-like  modes  as  TM  modes.  The  electric  field 
of  the  first  mode  in  the  source-free  region,  i.e.,  outside  the 
surface  S'(r),  is  determined  by  the  imhomogeneous  vector 
Helmholtz  equation  with  the  right-hand  side  expressing  the 
boundary  condition  on  S' 

(V2  +  /con2(r))E(r)  =  -tw/*oJ"(r)5(S'(r))  (1) 

where  w,  /to,  fc0  =  2tt/A,  and  n  are  the  angular  frequency,  the 
permeability  of  vacuum,  the  wavenumber  of  light  in  vacuum, 
and  the  refractive  index,  respectively.  The  equivalent  surface 


P 


Fig  1  A  schematic  representation  of  the  equivalent  diffraction  problem  as 
described  by  Love’s  principle,  (2)  and  (3).  The  surface  element  with  the 
unit  normal  n'  generates  the  optical  field  at  an  observation  point  P  in  the 
source-free  region,  i.e.,  outside  the  volume  V.  Helds  generated  by  a  real 
source  J  and  by  the  entire  surface  S'  are  identical.  The  point  labeled  O  is 
the  coordinate  system  origin. 

current  density  is  given  by  Love’s  expression 

J*(r)  =  2n'  x  H'(r)|s'  (2) 

where  H'(r)  is  the  magnetic  field  on  S'  produced  by  a  real 
source.  The  Dirac-delta  function  in  (1)  is  used  to  express  the 
electric  field  in  terms  of  Green’s  function  G.  The  volume 
integration  reduces  to  integration  over  the  surface,  where 
the  equivalent  surface  current  density,  (2),  is  known.  Con¬ 
sequently,  we  have  the  following  2-D  integral  representation 
for  the  TE-mode  fields  outside  the  surface  S' 

E(r)  =  2iujfio  f  dS'G( r,  r')n'  x  H^f')  (3) 

Js' 

and  H(r)  =  (l/iw/to)V  x  E(r).  This  equivalent  diffraction 
problem,  as  opposite  to  the  initial  problem  with  the  bulk 
source,  is  depicted  in  Fig.  1. 

The  fields  of  a  TM  mode  can  be  described  in  a  very  similar 
manner  by  using  the  appropriate  equivalent  surface  current, 
determined  by  Love’s  principle.  In  considering  a  TM  mode, 
we  use  a  bar,  for  example,  n  — >  n(r),  to  describe  possible 
birefringence;  e(r)  =  5 of*2 ( r),  E,  H,  J"  (for  TE  modes) 
_>  E,  H,  J*  (for  TM  modes).  The  equivalent  TM  surface 
current  becomes 

Js(r)  =  — 2n'  x  E'(r)|.S'  (4) 

and  the  TM-mode  fields  in  the  source-free  region 

H(r)  =  -2 iwe0J  dS'G( r,  r')n2(rV  *  E'(r')  (5) 

and  E(r)  =  (i/we0n2(r))V  x  H(r),  and  Green’s  function 
G(r,  r')  is  calculated  with  n(r). 

The  general  expressions  presented  above  for  the  fields  gen¬ 
erated  by  equivalent  surface  sources  can  be  simplified  in  case 
of  planar  waveguide  geometry  by  using  any  approximation 
allowing  for  the  separation  of  variables.  This  is  done  below  in 
two  steps:  first,  by  separating  basic  eigenfunctions  describing 
confinement  of  the  optical  field  in  a  planar  waveguide,  and 
then  by  calculating  the  residual  transverse  eigenfunctions  that 
match  the  actual  index  profile  in  the  plane  (see  further  discus¬ 
sion  in  Section  IV).  The  eigenfunctions  describing  confinement 
in  a  planar  waveguide  (z  is  the  vertical  coordinate)  are 
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determined  by  well-known  second-order  differential  equations 
with  the  boundary  conditions  for  the  waveguide-substrate  and 
the  waveguide-cladding  interfaces  (see  [1],  [4]),  with  (30  as 
the  propagation  constant  in  the  longitudinal  direction  (rc)  of 
a  slab  region.  The  solution  gives  eigenfunctions  and 

eigenvalues  of  the  propagation  constant  (3qu  with  the  mode 
number  v  equal  to  zero  or  a  positive  integer.  Green’s  function 
in  (3)  can  be  then  reduced  to  a  2-D  form 

G(r,  r')  =  Y,  G *(P,  *  =  0,  1,  2,  •  - 

V 

(6) 

(V^  +  ^)C?,(p^p,)  =  -^2(P”P,)  (7) 

where  V,,  is  a  2-D  gradient  operator  and  p  and  p'  are  the  2-D 
radius  vectors  in  the  plane. 

After  separating  the  longitudinal  and  vertical  coordinates, 
the  electric  and  magnetic  fields  at  the  input  surface  are 
determined  by  the  TE  eigenmodes  of  the  input  waveguide 

E'(r)|s'  =e,Eo<p5'(/)t>(c).  e'  =  z  x  n'  (8) 

H'(r)|.v  =  —  <PS'(l)  ( A'i  —7—^  +  (9) 

uj/lq  \  dz  J 

where  the  function  *  (/)  describes  the  field  distribution  on  the 
principal  surface  S',  l  is  a  transverse  coordinate  (see  Section 
IV),  the  index  v  and  a  possible  transverse  mode  quantization 
index  are  omitted  for  brevity,  z  is  a  unit  vector  along  the 
vertical  direction  and  e'  is  the  direction  of  the  input  electric 
field  vector,  and  the  field  amplitude  Eo  is  an  arbitrary  constant. 
To  calculate  the  fields  outside  the  surface  S'  according  to  (3), 
one  needs  the  explicit  form  of  the  equivalent  surface  current 
density 

J*(r)  =  2n'  x  H'(r)|.v  =  -e'  —  /w.v (/)</'(-)•  (10) 

o;//u 

The  general  surface-integral  representation  (3)  for  the  TE 
modes  can  be  thus  reduced  to  a  linear  integral  over  the 
transverse  coordinate  on  the  surface  S' 

E(r)  =  —  2i/)0 £<)'</' (z)  j  dl'e'Glp,  (II) 

where  lf  is  the  arc  length  along  the  input  boundary  in  the 
x-l  plane.  This  equation  describes  diffraction  between  the 
input  and  output  principal  surfaces  in  a  given  primitive  block 
(Fig.  2).  The  unit  vector  e'  appears  under  the  surface  integral 
as  the  consequence  of  the  vectorial  aspect  of  the  problem. 
Equation  (11)  is  used  in  the  next  section  to  calculate  the 
block  transfer  amplitude.  Representation  of  the  field  by  a 
linear  integral  is  convenient  for  planar  waveguide  problems. 
The  accuracy  of  this  representation  is  limited  by  the  degree 
of  vertical  confinement.  It  is  also  assumed  that  multiple 
reflections  of  the  diffracted  optical  field  from  the  output  and 
input  principal  surfaces  can  be  neglected.  This  is  the  case  with 
star  couplers  and  routers  (see  Sections  V  and  VI). 

There  is  a  point  to  be  emphasized  regarding  approximations 
being  used  in  this  approach.  The  use  of  Fresnel’s  electric 
and  magnetic  vector  decomposition  at  the  surface  and  Love’s 


S'  s 


<b) 


Fig.  2.  Schematic  representation  of  (a)  an  individual  primitive  block  showing 
the  general  geometry  of  a  star  coupler  and  (b)  an  integrated  optical  structure 
with  its  input,  output,  and  intermediate  principal  surfaces  separating  primitive 
blocks,  representative,  for  instance,  of  a  Dragone's  router.  X, j  stands  for  a 
number  of  waveguides  between  the  /th  and  ,/th  principal  surfaces. 


principle  for  the  equivalent  surface  sources  are  not  subject 
to  any  approximations.  On  the  other  hand,  while  the  sepa¬ 
ration  of  variables  is  an  approximation,  it  does  not  lead  to 
loss  of  accuracy  or  abandonment  of  the  vectorial  formalism. 
This  is  just  a  simplification  allowing  for  a  one-dimensional 
(1-D)  integral  representation  of  the  fields,  computationally 
convenient  for  planar  waveguides  with  good  confinement.  The 
vectorial  representation,  in  a  sense  of  being  able  to  calculate 
field  intensities  for  any  polarization  angle,  is  maintained. 

The  normalization  condition,  which  is  necessary  for  con¬ 
struction  of  a  transfer  amplitude,  requires  knowledge  of  the 
input  Poynting  vector  on  the  surface  S'.  According  to  (8)  and 
(9),  it  can  be  represented  as 

P'  =  ^  R<;(E'  x  H  *) 

=  ^-|^o|2ks'(/)|2 

2  cofiQ 

■  Refn'/JolVH^)!2  +  (12> 


or,  after  integration  over  the  surface  S', 

T>/  _  I  IP  |2 

~  2w//(l  1  ■ 


Here  the  asterisk  stands  for  complex  conjugation.  The  polar¬ 
ization  vectors  e  and  e '  are  assumed  to  be  real  and  describe 
linear  polarizations. 

In  the  case  of  TM  modes,  the  boundary  conditions  are 
different  and  there  is  possible  birefringence,  so  that  the  system 
of  equations  for  the  TM  eigenfunctions  of  vertical  quantization 
needs  to  be  modified  [3].  There  is  then  a  complete  one-to- 
one  correspondence  between  the  TE  and  TM  modes.  Only 
an  exchange  of  electric  and  magnetic  fields  and  the  use  of 
quantities  with  a  bar  are  required.  Therefore,  for  a  TM  mode, 
(8)-(10)  are  replaced  with 


H'(r)|S'  =S'Ho<Ps'W(z) 


(14) 
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E'(r)|s<  =  - 


H°  Vs'V) 


uj£0  n2(r) 


+  (15) 


dz 

J*(r)  =  — 2n'  x  E'(r)|s‘ 


2^°,  t  ffovMOV'M 

ue0n2(  r) 


(16) 


where  ()q  is  the  TM  mode  propagation  constant,  and  the  gen¬ 
eral  surface-integral  representation  (5)  for  the  TM  modes  can 
be  reduced  to  a  linear  integral  over  the  transverse  coordinate 
on  the  surface  S' 

H(r)  =  2iTioHoJ>(z)  f  dl' e'G(p,  p')<fiS'(1')-  (17) 

Js' 

The  magnetic  fieldjimplitude  H0  is  an  arbitrary  constant,  and 
Green’s  function  G(p,  p')  is  defined  with  equations  similar  to 
(6)  and  (7)  for  the  case  of  TM  modes. 

There  is  an  obvious  symmetry  (up  to  the  general  sign) 
between  the  propagation  properties  of  the  TE  and  TM  modes 
as  described  by  (II)  and  (17)  for  diffracted  fields.  In  the 
remainder,  unless  otherwise  indicated,  we  will  discuss  only 
the  TE  modes. 


III.  TRANSFER  AMPLITUDES 


These  expressions  for  outgoing  fields  on  the  surface  5  are 
very  similar  to  expressions  for  the  incoming  fields  on  the 
surface  S'  [see  (8)  and  (9)],  except  for  the  presence  of  the 
amplitude  transfer  factor.  Multiple  reflections  of  the  diffracted 
optical  field  from  the  output  and  input  principal  surfaces  are 
neglected  in  (18)  and  (19). 

The  Poynting  vector  of  the  outgoing  field  can  be  easily 
obtained  on  the  basis  of  (20)  and  (21) 

P  =  i  Ro.(E  x  H* ) 

=  \E0\2\<ps(1)\2\I(S,  S')|2 

2 w/<o 

•  R.c  ^n/)0|V>(z)|2  +  z*V>(z)  ( 2 2) 

or,  after  integration  over  the  surface  S, 

P««*  =  n^-\E0\2\t(S,S')\2.  (23) 

2t j//-o 

Thus,  the  quantity  \t(S.  S')|2,  where  l  is  defined  by  (18), 
represents  a  power  (or  intensity)  transfer  function  for  a  given 
primitive  block  (Fig.  1(a)] 

t(s.  s')  = n;  pV""1-  =  ms,  s>)\2-  <24> 

^  *  *totaI 


In  this  section,  we  derive  a  general  expression  for  the 
transfer  amplitude  based  on  (1 1).  The  field  E(r)  is  generated  in 
a  slab  region  by  an  effective  surface  source  on  S'.  It  propagates 
from  an  input  surface  5'  to  an  output  surface  S  according  to 
Huygens’  principle  and  enters  the  system  of  output  waveguides 
[see  Fig.  2(a)].  Let  us  suppose  the  total  eigenfunction  of  the 
output  waveguide  system  contains  a  product  of  a  new  unit 
vector  e  and  an  eigenfunction  describing  vertical  quantization 
This  representation  separates  the  vector  properties  and 
the  vertical  coordinate  dependence.  The  remaining  spatial 
dependence,  related  to  the  transverse  direction,  is  described 
by  an  eigenfunction  ips(l)  (Section  IV).  The  electric  field 
at  the  output  surface  is  given  by  a  product  of  the  vertical 
and  transverse  eigenfunctions  with  the  transfer  amplitude  that 
describes  propagation  and  diffraction  in  the  primitive  block. 

A  natural  way  to  define  a  dimensionless  transfer  amplitude 
for  a  primitive  block  bound  by  surfaces  S  and  S'  is  to  consider 
the  functional  and  vectorial  projections  of  the  electric  field 
E(r),  as  defined  by  (11),  onto  the  product  of  eigenfunctions 
of  the  output  waveguides,  divided  by  a  field  amplitude  Eq 

t(S ,  S')=YjsdSip*sm*(z)e  E(r)  (18) 

=  -2*ft>  [  dl  [  dZ'e-eVK0G(p,P>S'(i')-(19) 
Js  Js‘ 

Orthonormality  properties  of  the  eigenfunctions  ^  and  <p  were 
used  to  get  the  second  line  in  this  equation.  With  (18),  the 
electric  and  magnetic  fields  related  to  the  output  surface  S 
can  be  represented  as  follows: 

E(r)|s  =et(S,  S')Eo<ps(W>(z)  (2°) 

H(r)|s  =  t{S,  S')  ^  VsV)  (&i  ^  •  (21) 


The  total  transfer  amplitude  for  the  entire  optical  system 
[see  Fig.  2(b)]  is  represented  by  a  convolution  of  multiple 
transfer  amplitudes  of  primitive  blocks,  with  summation  over 
all  intermediate  principal  surfaces 

^in) 

=  ^  St, )  •••  S\)t(S\.  Sin)-  (25) 

s, .  Sj.  -  5„ 

The  meaning  of  the  sum  is  specified  below  (Section  VI).  The 
corresponding  transfer  function  for  the  entire  optical  system 
is  defined  as 

T(Simt.  Si,,)  —  |/-(Sou„  Si„)|  .  (26) 

At  this  point,  it  is  convenient  to  introduce  a  normalized 
transverse  eigenfunction  of  the  TE  diffracted  field,  according 
to  (19) 

<, p{l )  =  -2i/?0  I  dl' e  *  e'G(p,  p)vs’{l')-  (^7) 

It  is  easy  to  see  that  if  the  surfaces  S  and  S'  are  smooth 
and  have  low  curvatures,  which  is  the  case  of  real  integrated 
optical  circuits,  the  normalization  condition  for  the  function 
<p(l)  can  be  written  as 

[dim)  |2  =  Jsidl'\<pS'{r)\2  =  l -  (28) 

The  concept  of  a  normalized  transverse  eigenfunction  of  the 
diffracted  field  is  very  useful  because:  1)  the  function  <p(l) 
is  defined  uniquely  by  (27)  and  (28)  everywhere;  2)  it  can 
be  easily  generalized  to  the  case  of  several  primitive  blocks 
(Section  VI);  and  3)  the  transfer  amplitude  (19)  can  be 
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represented  as  a  complex  overlap  integral  of  two  normalized 
eigenfunctions  <p(l)  and  (ps(l)  on  the  principal  surface  S,  i.e., 

t(S.  S')  =  J^dl  <Ps(l)(fi(l).  (29) 

These  two  eigenfunctions  describe  the  incoming  and  outgoing 
fields  on  the  principal  surface  S . 

In  calculating  the  properties  of  optical  devices,  it  is  easy 
to  see  that  the  eigenfunctions  <ps(/)  and  ip(l)  should  closely 
match  each  other  on  every  principal,  or  output  surface  S', 
in  order  to  maximize  the  efficiency  [see  (24)  or  (26)]  of  a 
particular  device  design.  In  other  words,  each  pair  of  these 
two  functions  (for  a  given  S)  has  to  be  commensurable.  For 
example,  both  can  be  chosen  as  transverse  Bloch  functions  of 
the  same  order,  as  discussed  in  Section  VI. 

To  complete  our  consideration  of  transfer  amplitudes,  the 
corresponding  basic  equations  for  the  TM-mode  transfer  am¬ 
plitude  t(S ,  S"),  the  transverse  eigenfunction  y?(/),  and  the 
transfer  function  T(S.  S'),  are  presented  as 


7(.S\  S')  =  m,  /  <"  /  <U'e  e'vZ(l)C(p.  P'Ws.(l') 

Js  ./S' 

II 

•Gl 

■•A* 

% 

(30) 

<P(I)  -  1  <U'  e  •  eG(p,  p)ips,(l') 

(31) 

T(S,S')=  =|7(S,  S') |2. 

n  *  ‘■■total 

(32) 

All  the  quantities  in  (30)— (32)  are  marked  with  a  bar.  Because 
of  symmetry  between  the  two  sets  of  equations  (19),  (24), 
(27),  and  (29)  for  the  TE  modes  and  (30)-(32)  for  the  TM 
modes,  one  can  expect,  for  small  birefringence,  that  the 
transverse  eigenfunctions  and  transfer  amplitudes  for  the  two 
polarizations  will  be  numerically  very  close  (see  Section  IV). 

The  expressions  for  transfer  amplitudes,  transverse  eigen¬ 
functions,  and  the  block  and  total  transfer  functions  obtained 
in  this  section  represent  a  direct  consequence  of  vector  elec¬ 
trodynamics,  i.e.,  they  depend  on  the  light-polarization  vectors 
and  geometry  of  the  optical  system. 

IV.  Transverse  Eigenfunctions  for  Single 
and  Multiple  Waveguide  Structures 

The  transverse  coordinate  dependence  of  the  optical  field  on 
any  principal  surface  5,  described  above  by  an  eigenfunction 
ips{ /),  is  generally  not  known.  This  function  is  a  solution  of 
a  second  order  differential  equation  and  it  can  be,  in  prin¬ 
ciple,  calculated  exactly.  It  can  also  be  modeled  for  specific 
waveguide  designs  by  assuming  a  cosine,  or  more  complex, 
periodic  refractive  index  or  aperture  profiles  [14].  Any  reliable 
eigenmode  solver  can  be  used  with  our  proposed  method 
of  transfer  amplitudes.  We  describe  below  an  approach  for 
calculating  transverse  eigenmodes  that  is  similar  to  the  well- 
known  effective  index  method  in  respect  to  basic  equations, 
but  it  uses  a  different  numerical  technique  and  produces  some 
new  results. 

The  first  part  of  this  section  describes  a  simple  numerical 
method  for  calculating  the  transverse  eigenfunction  v?s(/).  I* 


is  assumed  that  the  vertical  quantization  is  performed  first 
and  that  the  index  profile  in  the  transverse  direction  is  an 
arbitrary  function  of  the  transverse  coordinate  (/).  We  are 
thus  able  to  use  vertical  eigenfunctions  to  obtain  ips{l)-  The 
propagation  constant  fto  is  interpreted  as  a  bare  propagation 
constant  corresponding  only  to  vertical  quantization  in  a  slab 
region.  In  case  of  an  arbitrary  index  profile  in  the  trans¬ 
verse  direction,  the  electromagnetic  field  dependence  on  the 
longitudinal  coordinate  (z)  is  determined  by  a  new  effective 
propagation  constant  ft,  to  be  determined  below. 

For  a  waveguide  system  of  arbitrary  geometry  in  the  x-l 
plane  (the  local  Cartesian  coordinates  are  x ,  /,  z)  with  strong 
confinement  in  the  z  direction,  a  general  TE  eigenmode  is  a 
solution  of  the  Helmholtz  equation  (V2  4-  A*oW2( r))E(r)  =  0, 
and,  far  enough  from  cutoff,  it  can  be  described  as 

E(r)  =  eE0eiflr<p(l)rl>(z)  ~  j**,  ft  £  (i0  (33) 

where  ft  and  ft o  are  the  effective~and  bare  propagation  con¬ 
stants,  respectively  (fto  enters  as  an  implicit  parameter  in 
V»(z)  and  </?(/);  fto  >  ft),  e  is  a  unit  vector,  and  /50  is  the 
initial  amplitude.  The  coordinate  /  reduces  to  the  Cartesian 
coordinate  y  when  the  principal  surface  S  is  a  planar  surface 
perpendicular  to  the  x  direction.  In  the  waveguide  plane,  the 
refractive  index  is  a  function  of  the  transverse  coordinate  7/.tr(/) 
with  the  maximum  value  n\ . 

The  difference  between  fto  and  ft  can  be  determined  as 
follows.  We  look  for  a  transverse  eigenfunction  ip(l)  describing 
the  waveguide  region  between  the  principal  surfaces,  so  that 
the  surface  index  ( S  or  S')  is  inappropriate  here.  After 
separating  the  longitudinal  coordinate,  the  product  of  rj)  and 
ip  in  (33)  should  satisfy 

{w>  +  lP>  +  (*g»2(Z.  z)  -  =  0  (34) 

where  u(L  z)  is  the  refractive  index  profile  in  the  transverse 
and  vertical  directions,  equal  to  nlr(/)  inside  the  layer.  In 
addition,  </;  must  satisfy  the  condition  of  vertical  quantization 

+  ^n2(z)  -  /ty’(z)  =  ()  (35) 

where  n(z)  is  the  “bare”  refractive  index  profile  in  the 
transverse  and  vertical  directions,  equal  to  n\  inside  the  layer. 
After  combining  (34)  and  (35),  the  basic  Schrodinger-like 
equation  for  the  eigenfunction  </>(/)  is  obtained 

(J^  +  A[e-V(l)\\V(l)=Q  (36) 

where  the  dimensionless  eigenvalue  5  and  the  ‘'potential  en¬ 
ergy”  of  a  photon  V(l),  which  are  limited  to  0  <  £,  V(l)  <  1, 
are  defined  as  follows: 

£  =  k2{n2  _  n2r(Z))  =  AV{1)'  (37) 

In  case  of  a  buried  waveguide(s),  when  the  transverse  index 
profile  ntr(l)  lies  in  between  n\  and  712,  the  parameter  A  = 
kl(n\  —  n\).  If  e  is  known,  then  the  unknown  ft  =  neff^o  = 
\/ fto  -  Ae  is  a  function  of  fto ,  A,  d,  a,  and  the  parameters  in 
V(l);  neff  is  the  effective  refractive  index.  In  this  approach. 
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Fig.  3.  Intensities  of  the  TE  and  TM  eigenmodes  of  a  single  waveguide  with 
a  rectangular  index  profile,  and  two  normalized  Bloch  functions  in  the  periodic 
waveguide  profile,  calculated  for  the  values  of  the  transverse  quasi-momentum 
q  =  0  and 

the  condition  of  vertical  quantization  is  included  explicitly 
through  (35). 

The  numerical  technique  based  on  (35)  and  (36)  is  ap¬ 
plicable  with  good  accuracy  to  any  index  profile  in  the 
transverse  direction,  including  profiles  with  abrupt  boundaries 
in  both  directions  and  a  smooth  profile  only  in  the  transverse 
direction.  The  procedure  of  vertical  quantization  within  the 
present  approach  is  exactly  the  same  as  in  the  effective  index 
method.  Some  examples  of  simultaneous  evaluation  of  the 
eigenfunction  y?(/)  and  the  eigenvalue  e  are  given  below. 

In  the  case  of  a  symmetric  slab  waveguide,  we  start  with 
the  known  transcendental  equation  for  flu  expressed  through 
a  dimensionless  variable  X 

tan(X)  =  1.  X  =  \  y(fcom)2  -  Jl  (38) 

where  U  =  {hk0/ 2)v/uf  -  uf  is  a  characteristic  parameter. 
This  equation  determines  fl o  as  a  function  of  A,  ni,  112 , 
the  layer  thickness,  and  a  quantum  number  of  the  verti¬ 
cal  quantization  v  —  0.  1.  2,  ••*.  In  the  case  of  a  single 
waveguide,  we  use  (36)  and  (37)  with  the  known  flo  to 
determine  eigenvalues  of  the  effective  propagation  constant 
1 3  and  the  corresponding  eigenfunctions  </?(/).  The  examples 
of  the  resulting  eigenfunctions  <p(l)  for  a  single  waveguide 
within  a  periodic  structure  with  the  refractive  index  difference 
An  =  0.0085  and  birefringence  0.0125  are  plotted  in  Fig.  3. 
A  rectangular  index  profile  is  assumed.  Both  the  TE  and  TM 
modes  are  shown.  The  Bloch  functions,  also  shown  in  Fig.  3, 
will  be  discussed  below.  There  is  a  strong  similarity  between 
different  eigenfunctions  in  Fig.  3,  especially  in  the  core  wave¬ 
guide  region.  For  example,  the  TE  and  TM  eigenmodes  almost 
coincide.  Our  numerical  calculations  have  also  shown  that  for 
any  refractive  index  profile  and  almost  everywhere  in  l  the 
transverse  eigenfunction  is  very  close  to  a  Gaussian  (see  also 
[11]). 

The  even  and  odd  TE  eigenmodes  of  a  symmetric  double 
waveguide  structure  in  the  case  of  a  strong  coupling  are 
plotted  in  Fig.  4.  The  presence  of  two  eigenfunctions  with 
different  parities,  and  having  slightly  different  eigenvalues,  is 


Fig.  4.  The  even  and  odd  TE  eigenmodes  of  a  symmetric  double  waveguide 
structure  with  strong  coupling. 


the  consequence  of  tunneling  through  the  central  barrier.  The 
tunneling  is  important  in  the  analysis  of  optical  devices  such 
as  3-dB  couplers  [1]. 

To  demonstrate  the  capability  of  our  numerical  method, 
we  use  (36)  to  evaluate  the  transverse  eigenfunction  of  a 
periodic  waveguide  structure  with  the  index  profile  given  by 
V(l)  =  c:os2(7r//r/),  where  d  is  the  period.  In  this  case,  the 
function  q>(l)  is  represented  by  a  periodic  Bloch  function 
which  can  be  expressed  through  standard  Mathieu  functions 
(see  [20,  ch.  20]) 


v>(0 


f pq(l)  =  cer 


(39) 


where  q  is  the  quasi-momentum  in  the  transverse  direction 
and  r  is  the  dimensionless  index  of  the  Mathieu  function. 
The  corresponding  eigenvalue  in  (36)  and  the  value  of  the 
propagation  constant  are  determined  as 


where  ar(Q)  is  the  characteristic  value  of  the  Mathieu  equa¬ 
tion  ce"(z,  Q)  +  (ar  -  2 Q  cos(2 z))  cer(z.  Q)  =  0.  For  the 
first  Brillouin  zone,  we  have  0  <  \q\  <  it/d  (0  <  \r\  <  1). 

The  normalizable  part  uq(l)  of  the  Bloch  function,  which 
according  to  Floquet’s  theorem  is  defined  by  the  equation 
c pq(l )  =  eiqluq(l ),  represents  an  element  of  the  orthonormal  set 
of  functions  f^<d^2  dl  u*(l)uq' (l)  =  8qq'.  This  normalizable 
part  is  plotted  in  Fig.  3  for  two  values  of  q  (0  and  7r/d). 
Within  the  waveguide,  all  solutions  are  very  close  to  each 
other.  For  q  =  7r/d,  the  Bloch  function  behaves  differently. 
The  periodic  Bloch  eigenmodes  obtained  above  cannot  be,  in 
general,  excited.  However,  they  can  be  excited  with  specific 
envelope  functions  and  specific  values  of  the  quasi-momentum 
q  (see  Section  V). 


V.  Fourier  Optics  Approximation  and  Star  Couplers 

In  this  section,  we  consider  propagation  of  eigenmodes  and 
transverse  eigenfunctions  based  on  the  idea  of  the  transfer 
amplitude  defined  in  (29).  The  general  formulation  discussed 
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in  Section  III  is  applied  here  to  efficiency  calculations  of:  1) 
Dragone’s  star  coupler  and  2)  waveguide  grating  router  (see 
Section  VI). 

It  is  well  known  that  diffraction  of  guided  waves  in  a  star 
coupler  with  circular  (cylindrical)  principal  surfaces,  periodic 
input,  and  output  waveguide  structures,  and  the  ratio  d/R  <C  1 
(where  d  is  the  distance  between  waveguides  in  the  transverse 
direction  and  R  is  the  characteristic  radius  of  principal  surfaces 
bounding  a  slab  region)  can  be  described  with  high  accuracy 
within  the  Fourier  optics  approximation  [14],  [211,  [22].  Cal¬ 
culations  of  transverse  eigenfunction  ip(l)  of  the  diffracted 
field  and  the  transfer  amplitude  t(S,  S')  for  a  star  coupler 
in  the  Fourier  optics  approximation  are  considered  below  for 
different  types  of  input  and  output  apertures.  Numerical  results 
for  the  transfer  function  are  presented  as  a  function  of  the 
total  number  of  waveguides  and  the  characteristic  aperture 
size.  However,  it  is  interesting  to  point  out  that  the  general 
method  described  here  is  applicable  to  non-Fourier  optics 
situations  as  well. 

The  explicit  form  of  the  2-D  Green's  function,  defined  by 
(7),  for  a  slab  region  of  a  star  coupler  [see  Fig.  2(a)]  is  given 
by  the  Hankel  function  of  the  first  kind  and  zeroth-order  (the 
quantization  index  in  the  propagation  constant  is  omitted  here 
for  simplicity) 


where  //^(.r)  =  •/(>(.'•*)  4-  /‘V„(t),  and  the  approximation  cor¬ 
responds  to  the  asymptotic  expansion  of  the  Hankel  function 
at  ,%\p  -  ff\  »  1  [ 20 1 .  The  distance  between  points  l}  and 
l}\  with  the  angular  coordinates  <\  and  n'  on  the  arcs  S  and 
S\  respectively,  is  determined  as  \p  -  p'\  =  l\{(r os  n  4- 
cos  <*'  -  I )’-  4-  (sin  n  -  sin  a/)J),/-\  where  the  characteristic 
values  of  the  input  and  output  aperture  angles  |n|.  |n'|  1 

[Fig.  2(a)].  Thus,  in  the  standard  diffraction  approximation 
[191,  we  substitute  |p  -  p'\  =  R  in  the  denominator  and 
\p  -  p'\  =  R :  (1  -  wv')  in  the  exponent  of  (41).  This  gives 
rise  to  the  propagation  function 


2  V  2n/i{)R 


./Til  /?( 1  — 


where  we  explicitly  use  the  angular  variables  o  and  o'  instead 
of  a  pair  of  2-D  radius  vectors  p  and  p' .  It  is  also  convenient 
to  normalize  the  eigenfunctions  by  the  angle  o  instead  of  the 
arc  length  l  =  aR.  [see  (28)] 


J  <li\<p( Ol2  =  J  d<A<Kn)\2  =  i-  ¥>(0  =  -J=i  (K")- 

(43) 

For  a  TE-polarized  eigenmode,  the  product  of  the  initial  and 
final  polarization  vectors  reduces  to  unity  for  |a|,  |a'|  <  1 
[see  (19)  and  (27)]:  e  •  e'  =  cos  a  cos  a'  —  sin  a  sin  <v'  ^  1. 
The  basic  equations  (27)  and  (29)  can  be  reformulated  as 

t(S,S')  =  Jda<fs(a)<j>(a)  (44) 


where  the  integral  in  (45)  is  a  typical  Fourier  optics  integral. 
For  diffraction  problems  under  consideration,  numerical  cal¬ 
culations  show  that  the  difference  between  the  results  obtained 
with  the  approximate  propagation  function  (42)  and  the  exact 
expression  (i/4)H^\fio\p  ~  p'\)  is  negligible. 

The  oscillating  nature  of  the  diffracted  field  on  the  output 
surface  5  is  important.  If  the  input  into  the  slab  region  of 
the  coupler  is  introduced  through  a  single  waveguide  with  the 
ordinal  number  in  on  S',  the  transverse  eigenfunction  (45)  of 
the  diffracted  field  on  S  can  be  represented  in  general  by  the 
Bloch  function  of  the  order  m  (m  =  0.  ±1.  ±2.  • ): 

=  (46) 

The  exponent  in  (46)  can  be  written- in  the  real  coordinate  form 
as  imftodfY  —  />//.  where  /  =  n/t  and  <j  —  mS{id/R  is  the 
transverse  Bloch  momentum.  Similar  Bloch-mode  excitations 
have  been  discussed  by  Dragone  in  the  case  of  merging  or 
tapering  multiwaveguide  structures  1 16].  The  angular  period  of 
Bloch  oscillations  equals  27r/(|m|/W)  and  is  much  smaller, 
for  |m|  >  1,  than  the  characteristic  spread  of  the  diffracted 
field  described  by  the  function  f/>„(n).  This  spread,  also  called 
the  diffraction  angle  o^h,  is  of  the  order  It r/(/V/),  as 
determined  below.  On  the  other  hand,  this  period  is  much 
greater  than  the  natural  angular  period  of  the  waveguiding 
system  d/R.  Consequently,  we  obtain  the  following  sequence 
of  inequalities  (|m|  >  1): 


,,J=«<'w  =  i4W<n . =  'u  {  n 

Diffracted  waves  generated  by  a  side- waveguide  input  in  a 
circular  star-coupler  geometry  can  be  approximated  by  Bloch 
modes  described  by  periodic  boundary  conditions.  In  this  case, 
we  have  [see  (39)  and  (40)  and  discussion  of  the  Bloch 
functions  in  Section  IV |:  <j  ~  27rm/Ay/  -  tii.i{)d/R  = 
(n/d)r,  i\Tt  —  2irR/(iudr  ~  10,  r  —  mS{)dn\i/~  —  — 

0.2m,  and  [m\  <  5  for  the  first  Brillouin  zone.  The  periodic 
boundary  is  applied  at  a  distance  Ntd. 

It  should  be  emphasized  that  in  the  overlap  integral  in  / 
[see  (44)1,  Bloch  oscillations  in  </>(n)  cannot  be  compensated 
by  a  choice  of  the  appropriate  value  of  d/R  because  m  is 
a  variable  (m  =  0.  ±  1 .  ±2.  • ).  For  the  central  waveguide 

input,  the  Bloch  factor  equals  !(///=  0).  The  field  produced 
on  the  output  surface  S  by  each  side-waveguide  input  (m  ^  0) 
is  harmonic  [see  (46)].  This  field  enters  the  output  waveguide 
structure  and  propagates  through  it,  conserving  the  specific 
phase  relationship  (for  the  transverse  Bloch  functions)  between 
different  waveguides  (see  Section  VI).  Consequently,  to  pro¬ 
vide  the  highest  overlap  of  the  transfer  amplitude  [see  (44)], 
the  eigenfunction  at  the  output  determined  on  5,  should 

contain  a  Bloch  function  factor  (e“,m'*°'/a)  as  well,  with  the 
same  value  of  rn  as  in  <fi(a)  [see  (45)  and  (46)]. 

For  a  single-waveguide  input  with  a  transverse  rectangular 
aperture  2a,  the  nonperiodic  part  <£0  of  the  eigenfunction  is 
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Fig.  5.  Intensities  of  the  nonperiodic  parts  <&)(«)  of  transverse  eigenfunc¬ 
tions  as  determined  from  (48)  (2 a  =  6  //m)  and  (49)  ( d  =s  10  /xm).  The  output 
intensity  generated  by  a  Gaussian  input  eigenmode,  as  an  approximation  to 
the  TE  and  TM  eigenmodes  in  Fig.  3,  is  also  shown. 

given  by  (45)  as 

U»)  =  (Z^)  l/*e<<,il«nc  (fiaa)  (48) 

where  sinc(:r)  =  sin(:/;)/:/:.  There  are  two  special  cases  of 
a  =  d/ 2,  which  corresponds  to  the  case  of  waveguides 
merging  at  S',  and  a  =  dj 4,  which  corresponds  to  a  periodic 
step-index  waveguide  structure  on  S'. 

For  a  single  waveguide  input  with  a  periodic  cos2(7 rRa/d) 
aperture  profile,  the  nonperiodic  part  <j) o  of  the  eigenfunction 
is  given  by 


The  first  zero  of  this  function  occurs  at  a  =  ocdiff  =  47r//3od. 
The  resulting  angular  distributions  of  along  with  the 

distribution  generated  by  a  Gaussian  input  eigenmode  are 
depicted  in  Fig.  5. 

When  the  eigenfunctions  and  <f>s{<*)  are  known, 

the  integral  in  (44)  can  be  easily  evaluated.  The  calculated 
transfer  functions  T  —  \t\2  are  presented  in  Fig.  6  for  different 
values  of  N.  Transfer  functions  are  sensitive  to  the  size,  the 
total  number,  and  the  periodicity  of  the  input  and  output 
waveguides.  This  allows  for  easy  optimization  of  the  device. 
The  calculated  values  for  the  star  coupler’s  efficiency  are 
realistic  and  relatively  high,  reaching  a  maximum  of  ~80%. 
This  efficiency  is  obtained  with  the  condition  of  zero  aperture 
outside  the  waveguide,  as  used  previously  by  Dragone  [14]. 

VI.  Simulation  of  a  Waveguide  Grating  Router 

Having  calculated  the  diffracted  field  distribution  in  a  star 
coupler,  we  now  turn  to  a  more  complicated  optical  system 
of  two  star  couplers  connected  by  a  multiwaveguide  struc¬ 
ture,  a  version  of  Dragone’ s  router  [15],  [23]  (Fig.  7).  The 
discussion  of  a  waveguide  grating  router  consists  of  three 
parts.  We  first  discuss  a  connection  between  the  convolution  of 
transfer  amplitudes  and  Bloch  modes.  A  zero  waveguide-width 
approximation  is  then  used  to  obtain  relationships  between 


Fig.  6.  Efficiency  T(S,  S')  of  a  star  coupler  plotted  as  a  function  of  the 
number  N  of  output  waveguides  in  the  case  of  a  single-waveguide  input.  We 
show  two  cases  of  a  waveguide  system  merging  on  5  with  a  ~  d/2  and 
a  =  d/4  [(48),  solid  and  dashed  lines,  respectively],  and  a  case  of  a  periodic 
cos2  ( x Ra/d)  aperture  profile  on  S  [(49),  dashed-dotted  line]. 


Fig.  7.  An  optical  system  consisting  of  two  star  couplers  connected  by  the 
multiwaveguide  structure,  known  as  the  Dragone’s  router. 

various  waveguide  parameters  of  the  router.  In  the  third  part, 
we  calculate  the  optical  response  function  of  the  router  and 
present  results  of  numerical  simulations. 

In  the  router  illustrated  in  Fig.  7,  there  are  four  principal 
surfaces:  the  input  surface  5in  with  the  angular  variable  a', 
the  two  intermediate  surfaces  Si  and  S[  (variables  ai  and  a\, 
respectively),  and  the  output  surface  Sont  with  the  variable 
a.  For  this  system  of  surfaces,  the  direct  application  of  (25) 
results  in  the  following  specific  form: 

t(Sout,  Sm)  =  S{)t{Sit  ^in)  (^0) 

s,.s; 

where  t{Su  Sin)  and  t{SouU  S[)  are  the  transfer  amplitudes 
for  the  two  star  couplers,  i.e„  the  two  slab  regions;  and 
t(S[,  Si)  is  the  transfer  amplitude  for  the  waveguide  grating. 
However,  as  written,  this  is  only  a  formal  relation  and  the 
meaning  of  summations  over  principal  surfaces  is  not  defined. 

It  is  convenient  to  use  equations  similar  to  (27)  and  (29) 
[or  (44)  and  (45)  for  the  angular  variables]  to  express  the  total 
transfer  amplitude  t(S0 ut,  S)n)  in  terms  of  diffracted  fields  in 
Dragone’s  router.  From  the  purely  geometrical  consideration, 
one  can  establish  that  the  transmission  through  the  multi- 
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waveguide  grating,  from  5X  to  S[  [represented  formally  by 
/.(S',  Si)  in  (50)],  can  be  described  by  a  specific  transfer 
amplitude,  called  here  W{ci\,  ai),  which  is  related  to  the 
grating  geometiy  and  the  detailed  shape  of  the  waveguide 
openings  at  the  two  surfaces.  When  both  star  couplers  are 
of  the  same  shape  and  size,  the  transfer  amplitude  W  is 
proportional  to  the  delta  function  of  its  arguments,  and  we 
get  the  following  symmetric  form: 

W((y\,  r*i)  =  -  «i)  (51) 

where  the  single-parameter  phase  factor  tr(n  i )  describes  phase 
shifts  between  adjacent  waveguides  arising  from  their  relative 
displacements  and  different  optical  paths  between  Si  and  SJ. 
Since  the  radiation  and  absorption  losses  are  assumed  to  be 
negligible,  the  absolute  value  of  the  phase  factor  is  taken  as 
1.  For  the  two  star  couplers  of  different  sizes,  a  scaling  factor 
should  be  introduced  in  (51). 

The  Dirac-delta  function  in  (51)  is  very  important  since  it 
can  be  easily  expanded  into  complex  Fourier  series  |24] 

-  fV|)  =  Al </A(/V/(o'i  -  €1  I  }) 


where  it  is  easy  to  recognize  the  Bloch  exponentials  similar  to 
that  defined  in  (46).  The  product  0()  d  is  simply  a  normalization 
factor.  As  a  result,  the  waveguide  transfer  amplitude  can  be 
decomposed  in  terms  of  a  complete  orthogonal  set  of  basic 
Bloch  modes. 

Consequently,  if  we  interpret  the  summation  in  (50)  as 
integration  over  the  angular  coordinates  on  each  principal 
surface,  that  equation  can  be  reformulated  as 

•"X; 

S-i„)=  Y  ■■)  (53) 

III  —  -  X 

which  means  that  the  Bloch  function  with  a  specific  number  m 
propagates  through  the  waveguide  grating  without  distortion. 
Partial  transfer  amplitudes,  for  each  Bloch  index  w,  are 
determined  by  substituting  the  formal  summation  over  S\  or 
,5’j  with  integration  over  rvi  or  n\.  This  representation,  with 
partial  transfer  amplitudes  for  each  index  ///,  (53),  contains 
only  double  diffraction  in  the  two  star  couplers.  The  Bloch¬ 
mode  transmission  through  the  grating  is  only  implied.  While 
in  this  shortened  description  there  is  no  explicit  summation 
over  surfaces  .S' j  and  $[,  this  picture  is  still  equivalent  to  the 
initial  (50). 

Once  the  summation  procedure  over  principal  surfaces  is 
specified,  (50)  can  be  reduced  to  the  integral  form  (44).  Using 
(51),  the  general  formulas  describing  double  diffraction  in  two 
star  couplers  can  be  represented  as 

*(£.„»,  sin)  =  jf  d«rs . («)*(«)  (54) 

<j>(rt)  =  —  2i(JoR  j  (ln\  G(n,  (x\ )  daiw(cti) 


x  f>{n\  -  «, )  (-2 :ipoR)  /  da'  G(a,.  a') 

Jsin 

•0slB(a')  (55 


or,  after  introducing  the  integral  kernel  A\ 

(f){a)  =  J  dafK(a,  a')0s, \a(nf)  (56) 

K(a,  af)  =  (-2ift0R)2  f  daiG(a,  ai)w((\x)G(ai.  a). 

d  Si 

(57) 

The  remaining  double  integral  in  (56)  and  (57)  is  quite 
simple  in  the  Fourier  optics  approximation.  Let's  consider 
a  zero-width  approximation  in  which  the  grating  produces 
a  periodic  signal  on  the  output  surface  The  charac¬ 

teristic  function  w  in  (55)  can  be  represented  by  a  sum  of 
discrete  contributions  from  each  waveguide  in  the  grating; 
each  contribution  is  included  with  a  phase  factor  describing 
the  corresponding  optical  path 

«»(fti)  =  <u  y  -  «»rf)  <5g) 

» 

where  A  is  the  length  difference  between  adjacent  waveguides 
connecting  S \  and  S\,  the  product  of  /7A  is  the  corre¬ 
sponding  optical  path  difference  between  waveguides,  and 
s  =  0.  ±1.  ±2,  •  *.  The  factor  n(/  was  chosen  here  for 
convenience.  By  using  (42)  for  G(n.  nf)  and  performing  the 
integration  over  infinite  limits,  we  get  an  analytical  expression 


K(n,  a')  =  -i(10[{c2ll1“n  ^  Y  ,:'S 


s(/JA-;V/(oh 


The  complex  Fourier  series  can  be  expressed  in  terms  of  a 
periodic  delta  function  using  the  well-known  formula  related 
to  Floquet’s  theorem 


^  v.  ',s^  =  27 r  ^  />(£  +  2 7T/)) 


where  ^  is  a  variable  and  .s,  p 
(59)  is  then  represented  as 

K{ a,  n)  =  -ic2iii"R  V  A 


(I.  ±1.  ±2. 


Equation 


fi  A  2t tj) 

o  +  o  -  —  —  + 

A)  d  A) <7 


The  grating  thus  generates  a  periodic  delta  function  in  the 
transfer  amplitude.  The  argument  of  the  delta  function  means 
that  the  output  field  in  (56)  is  a  set  of  shifted  and  inverted 
images  of  the  input  field.  In  practice,  the  integration  is  per¬ 
formed  over  a  finite  region  and  a  smooth  bell-shaped  function 
substitutes  for  each  delta  function  in  (61). 

Equations  (56)  and  (61)  have  important  applications  to 
wavelength  multiplexing  and  demultiplexing  in  integrated  op¬ 
tical  circuits.  Suppose  the  input  to  the  router  consists  of 
four  closely  spaced  and  equidistant  wavelengths  in  a  single 
waveguide,  with  four  corresponding  values  of  propagation 
constants  ft  and  ftQ.  This  set  of  wavelengths  can  be  represented 
by 


A n  —  Ao(l  T  ns) 

A.  =0(1  -  ne), 


A >»»  =  0o  (1  -  ne) 


where  7t  —  0,  1,  2,  3,  and  A0,  0,  and  ft0  are  the  initial  values 
for  the  (de)multiplexing  problem.  These  expressions  for  ft  and 
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p{)  are  valid  for  e  <  1.  To  resolve  two  adjacent  wavelengths, 
the  following  condition  on  the  last  term  in  the  argument  of  the 
delta  function  in  (61)  must  be  satisfied: 

2***  ,27 m=2KPoe  =  nd  m 

fihn+id  (hud  pod 

where  pa  is  the  order  of  the  diffraction  maximum.  This 
condition  includes  contributions  from  propagation  constants 
of  different  wavelengths.  Higher  order  maxima  are  present, 
with  the  same  intensity,  due  to  the  periodicity  of  the  delta 
function  in  (61).  Thus,  if  the  parameter  p()  is  chosen  in  such 
a  way  as  to  make  (64)  equal  to  the  angular  period  afi  of  the 
waveguide  structure,  then  the  o'  in  the  delta  function  in  (61) 
takes  the  following  values: 


the  efficiency  of  a  router:  the  efficiency  calculated  in  this  way 
is  always  underestimated.  However,  (64)  and  (66)  represent  a 
general  relationship  between  waveguide  parameters  and  are  in 
fact  approximation-independent. 

To  simulate  a  real  router  with  a  finite  waveguide  width, 
we  use  the  idea  of  universality  of  diffracted  eigenmodes  in 
star  couplers.  For  identical  star  couplers,  the  nonperiodic  parts 
of  diffracted  eigenmodes,  </>t,( a),  generated  by  single  opening 
apertures  are  exactly  the  same,  see  (46).  This  has  two  direct 
consequences.  First,  the  diffracted  wave  in  the  first  star  coupler 
gives  rise  to  the  interception  factor  fs  in  the  optical  response 
function  of  a  router.  The  factor  fs  is  defined  for  each  input 
wavelength  and  each  waveguide  in  the  grating  as 

1/2 


rS<\,l  +0(1 


n  =  — a 


na<i  -1- 


(65) 

/-( 

\Ai  <i  *  J 

’  dee'j 


(68) 


where  n  =  0.  1.  2.  3.  The  light  of  any  wavelength  can  be  thus 
directed  to  any  output  waveguide,  depending  only  on  the  input 
waveguide  used  and  the  parameter  A.  The  phase  shift  a<//c 
in  the  argument  of  (65)  can  be  partly  or  entirely  compensated 
for  by  choosing  the  value  of  the  characteristic  path  difference 
A.  For  example,  complete  compensation  is  achieved  for 


=  0.  ±1,  ±2,  *  •  The  factor  fs  is 
a  smooth  function  of  s  because  the  angular  width  of  each 
waveguide  o  j  is  much  smaller  than  the  diffraction  angle  «fnir 
in  (47).  The  interception  factor  is  a  measure  of  the  efficiency 
of  a  single  star  coupler 


7 


A  =  ^rZ. 

t>  e 


l-LI2 


(69) 


(66) 


The  light  of  wavelength  A()  entering  the  central  input  wave¬ 
guide  is  thus  directed  to  the  central  output  waveguide.  The 
same  light  entering  the  waveguide  m  is  directed  to  the  output 
waveguide  —///.  The  light  of  wavelength  Ai  is  displaced  by 
one  angular  period  of  the  grating,  etc. 

The  different  wavelength  contributions  are  totally  separated 
(demultiplexed)  in  space  (angular)  variable  and 


ri(n)  ~  ^  n  - 


(67) 


(see  also  the  last  paragraph  of  Section  V).  By  varying  the 
waveguide  width  (2a)  and  the  total  number  of  the  output 
waveguides  (Ar)  of  the  star  coupler,  the  efficiency  can  be 
optimized  (see  Fig.  6). 

Second,  according  to  (46),  the  nonperiodic  diffracted  ampli¬ 
tudes  <!>[)(< \ )  produced  by  each  input  opening  in  the  second  star 
coupler  are  also  identical.  Thus  the  router  response  function 
itself  is  proportional  to  </> o(n  ).  Collecting  these  ideas  together, 
we  can  construct  the  following  expression  for  the  normalized 
amplitude  response  function  of  a  router  with  a  finite  waveguide 
width: 


where  r„  is  the  amplitude  coefficient  of  a  contribution  with 
A„  in  the  input  signal.  We  conclude  that  the  input  signal 
is  totally  demultiplexed  on  the  output  surface  of  Dragone’s 
router,  i.e.,  each  wavelength  enters  its  own  output  waveguide. 
When  the  phase  shift  in  (65)  is  not  completely  compensated 
for,  an  additional  phase  shift  should  be  entered  into  (67).  This 
might  be  useful  when  it  is  desirable  to  send  any  wavelength 
to  any  particular  output  waveguide. 

It  is  easy  to  estimate  the  angular  dimension  of  the  grating 
by  using  real  values  for  the  initial  wavelengths  and  the 
channel  separation:  for  instance,  A0  —  1.552  52  //xn  and 
AA  =  0.0032  //m,  respectively.  In  this  case,  e  =  AA/A0  = 
2.06  x  HP3.  This  corresponds  to  the  diffraction  maximum 
of  the  order  p{)  =  P0  d.a(l/2nc  -  106  in  (67).  By  making 
the  optical  path  differences  A  equal  to  (pd\,  where  d\  is  the 
spatial  period  of  waveguides  in  the  grating,  and  compensating 
for  the  phase  difference  of  (66),  we  obtain  ip  —  a^d/edi.  For 
d}  =  100  //in.  ten  times  larger  than  the  waveguide  period  d  at 
the  principal  surface,  we  calculate  ip  —  27.75°. 

Equations  (64)  and  (66)  were  obtained  in  the  zero- width 
approximation.  This  approximation  cannot  be  used  to  calculate 


. .  (70) 

where  the  summation  over  .s  is  carried  out  for  all  waveguides  in 
the  grating  and  m  is  the  ordinal  number  of  the  input  waveguide 
of  the  router.  According  to  (64)  and  (66),  the  combination  p A 
in  the  exponent  equals  27rpo  where  the  order  of  the  diffraction 
maximum  po  is  integer. 

An  example  of  wavelength  demultiplexing  based  on  (70) 
is  shown  in  Fig.  8  which  plots  the  diffracted  intensities  for 
four  wavelengths  as  a  function  of  the  output  angle.  Numerical 
results  for  wavelength  multiplexing  produced  in  the  central 
channel  by  the  two  adjacent  waveguide  inputs  with  two 
different  wavelengths  are  illustrated  in  Fig.  9  for  the  total 
number  of  waveguides  in  the  grating  N  29  and  2a  —  6  //. m. 

A  simple  estimate  can  be  given  to  show  the  kind  of  effi¬ 
ciencies  expected  from  this  description  of  a  router.  Assuming 
a  rectangular  shape  of  the  input  eigenmode  <j>sXn  («'),  we  make 
use  of  the  sinc-function  representation  [see  (48)]  and  approxi¬ 
mate  the  interception  factor  by  using  the  mean-value  theorem 
for  the  integral  [see  (68))  fs  ~  yri^|</>n(.SY¥,/)|.  A  rough 
estimate  of  (70)  can  be  obtained  by  changing  the  summation 
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Fig.  8.  The  output  intensity  of  a  router  as  a  function  of  the  angular 
coordinate;  wavelength  demultiplexing  produced  by  a  multiwavelength  input 
at  the  central  waveguide. 


Fig.  9.  An  example  of  wavelength  multiplexing:  inputs  into  two  adjacent 
waveguides  result  in  the  output  at  the  central  waveguide.  The  output  intensity 
(on  a  decibel  scale)  is  plotted  as  a  function  of  wavelength.  The  input 
wavelengths  are  A»  into  the  central  waveguide  and  Af>  +  A  A  into  the  first 
adjacent  side  waveguide. 


to  integration  over  «s  and  using  the  Fourier  transform  of  a  sine 
function  et/txsiii( ■(:/;)  da  =  i r  0(1  —  /r),  where  0(:r)  is 

the  step  function).  The  result  can  be  represented  as 

(  4t ru  \  1/2 

</>(«)  ~  0(«)  (  J  “  ran,,).  (71) 

In  this  approximation,  the  output  signal  is  also  rectangular 
with  the  width  equal  to  2a.  Because  of  the  strong  inequality 
af/  <  the  peak  height  is  proportional  to  0o(-ma,/), 

i.e.,  it  slowly  decreases  with  |m|.  The  maximum  intercept  is 
provided  in  the  case  of  waveguides  merging  on  the  output 
surface  of  the  first  star  coupler  [a  =  d/ 2,  see  (69)].  Using  the 
normalization  condition  for  the  eigenmode  0siti,  the  estimated 
integrated  peak  intensity  of  the  output  signal  (71)  is  given  by 

J  \$(a)\2  da  =  jj^\<l>o(-Tnad)\2 .  (72) 

According  to  (48),  |0o(O)|2  is  proportional  to  a.  In  the  case 
of  a  central  input  and  output  (m  =  0),  we  get  the  following 
numerical  values  of  the  router  efficiency:  60%  for  a  =  3  /z m 
and  80%  for  a  =  4  /zmi,  both  for  d  —  10  ^m. 


In  reality,  the  shape  of  the  input  eigenmode  0s5n  (<*')  and  the 
router  amplitude  response  function  0(a)  are  not  rectangular 
and  there  are  no  analytical  expressions  for  0(a).  Numerical 
calculations  show,  however,  that  we  can  still  use  an  expression 
for  the  response  function  similar  to  that  given  by  (71).  In 
general,  the  effect  of  demultiplexing  is  described  by 

3 

0(a)  ~  <f>o{-'mad)Q  ^  cn0^(-a  -  (m  +  n)a</)  (73) 

n— 0 

with  the  obvious  features  of  (67)  and  (71).  Here  Q  is  the 
numerical  factor  which  reduces  to  the  square-root  factor  in 
(71)  in  the  rectangular  approximation,  and  <t>dd{a)  is  the 
exact  double-diffracted  transverse  eigenmode  on  the  output 
surface  of  the  router.  The  router  efficiency  can  be  obtained  by 
integration  of  \4>dd\2  over  the  angular  aperture  of  a  given  output 
waveguide.  In  this  way,  we  can  get  the  following  numerical 
values  of  the  router  efficiency:  54%  Jor  a  =  3  //m  and  72%  for 
a  =  4  /tm  which  are  in  reasonable  agreement  with  the  results 
of  the  rectangular  approximation  (2a  here  is  the  waveguide 
width). 

The  resolving  capability  of  the  device  depends  on  the  total 
number  of  waveguides  in  the  grating  section.  On  the  other 
hand,  the  total  efficiency  is  determined  by  the  waveguide 
width.  Fig.  10(a)  and  (b)  shows  two  examples  of  optimization 
of  a  waveguide  grating  router,  one  for  2a  =  G  // m  and  the 
efficiency  54%  and  another  for  2a  =  8  /im  and  the  efficiency 
72%,  as  a  function  of  N .  The  principal  maxima  in  both 
figures  change  only  slightly  while  the  first  sidelobes  decrease 
dramatically  with  increasing  N.  In  the  case  of  N  =  29 
and  2 a  =  8  //.m  [Fig.  10(b)l,  the  difference  between  the 
principal  maximum  and  the  first  sidelobe  exceeding  30  dB 
is  demonstrated  in  this  simulation. 


VII.  Summary 

We  present  a  method  of  modeling  and  simulating  inte¬ 
grated  optical  circuits.  It  is  based  on  the  vector-electrodynamic 
formulation  of  field  distributions  and  transfer  amplitudes  in 
optical  waveguides,  and  it  is  equivalent  to  Huygens’  picture  of 
propagation  and  diffraction  of  time-harmonic  electromagnetic 
fields.  Thus,  it  has  a  clear  physical  background.  The  method 
consists  of  the  following  steps: 

1)  dividing  a  waveguide  system  into  a  set  of  primitive 
blocks  separated  by  principal  optical  surfaces,  and  de¬ 
termining  boundary  values  of  electromagnetic  fields  and 
their  polarization; 

2)  expressing  the  diffracted  fields  through  surface  inte¬ 
grals,  generated  according  to  Love’s  field  equivalence 
principle,  by  the  equivalent  surface  current  densities; 

3)  computing  the  transfer  amplitude  for  each  block  sepa¬ 
rately,  and  accounting  accurately  for  effects  of  diffrac¬ 
tion; 

4)  calculating  the  transfer  amplitude  of  the  entire  optical 
system  as  a  convolution  integral  of  transfer  amplitudes 
of  its  primitive  blocks; 

5)  evaluating  the  total  transfer  amplitude,  as  a  function  of 
physical  characteristics  of  a  device,  such  as  the  size  and 
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(a) 


Fii:.  10.  Resolving  capability  of*  a  grating  router  lor  the  total  number  of 
waveguides  .Y  =  19.  23.  29  in  the  grating:  (a)  2a  =  ('*  //in  and  (b) 
2«i  —  «S  //  m. 


shape  of  apertures  on  principal  surfaces  or  the  number 
of  waveguides  in  the  grating. 

Fast  efficiency  optimization  of  each  primitive  block,  and 
the  entire  optical  system,  is  made  possible  in  this  method  by 
combining  analytical  and  numerical  calculations.  It  is  shown 
that  the  concept  of  a  normalized  transverse  eigenfunction  of 
the  diffracted  field  is  very  useful  in  calculations  of  transfer 
amplitudes.  The  crucial  role  of  transverse  Bloch  modes  in 
propagation  through  multiwaveguide  structures  is  emphasized. 
The  accuracy  of  this  method  is  demonstrated  by  calculat¬ 
ing  two  examples  in  detail.  First,  the  transverse  eigenvalues 
and  eigenfunctions  of  an  arbitrary  waveguide  structure  are 
described  by  the  basic  Schrodinger-like  equation,  combined 
with  the  standard  vertical  quantization  procedure.  The  general 
approach  is  also  illustrated  with  the  transfer  amplitude  and 
efficiency  calculations  for  a  Dragone’ s  star  coupler  and  router. 
The  method  is  applicable  to  a  variety  of  materials  and  systems. 

In  conclusion,  the  method  of  representing  complex  optical 
systems  developed  in  this  paper  provides  a  quantitative  de¬ 
scription  of  the  propagation  of  eigenmodes  and  is  a  powerful 
tool  for  numerical  calculations  of  diffracted  field  distributions, 
efficiencies,  etc.,  needed  for  complete  characterization  of  such 
systems. 
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We  report  on  low-frequency  noise  characteristics  of  visible-blind  GaN  p-n  junction  photodetectors. 

Carrier  hopping  through  defect  states  in  the  space  charge  region,  believed  to  be  associated  with 
dislocations,  is  identified  as  the  main  mechanism  responsible  for  the  dark  conductivity  of  the 
photodiodes.  Under  reverse  bias,  the  dark  current  noise  has  the  1//  character  and  obeys  the  Hooge 
relation  with  a~  3.  Under  forward  bias,  we  observe  generation-recombination  noise  related  to  a  trap 
level  with  the  activation  energy  of  0.49  cV.  Under  illumination,  detectivity  is  found  to  be  shot  noise 
limited.  The  noise  equivalent  power  of  a  200X200  /-an2  photodetector  is  estimated  at  6.6 
X 10“ 15  W/Hz,/2  at  a  bias  of  -3  V.  ©  1998  American  Institute  of  Physics . 
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I.  INTRODUCTION 

Ultraviolet-sensitive  photodetectors  with  a  sharp  cutoff 
for  the  visible  spectral  range  have  a  number  of  applications 
in  spectroscopy,  flame  sensing,  etc.  Due  to  the  large  direct 
band  gap,  gallium  nitride  is  ideally  suitable  for  the  prepara¬ 
tion  of  such  photodetectors.  Significant  progress  has  already 
been  achieved  in  the  development  of  Schottky  barrier  GaN- 
based  detectors,1,2  resulting  in  responsivities  as  high  as  0.18 
AAV  and  a  bandwidth  in  excess  of  10  MHz.  However,  the 
performance  of  Schottky  diodes  is  typically  limited  by  rela¬ 
tively  large  leakage  currents.  The  reverse  bias  leakage  cur- 
'  rent  in  a  p-n  junction  diode  is  expected  to  be  significantly 
lower  than  in  a  Schottky  barrier  diode  due  to  the  larger  bar¬ 
rier  height  This  promises  a  corresponding  performance 
advantage.3,4  Indeed,  p-n  junction  based  GaN  photodiodes 
with  dark  currents  as  low  as  a  few  pA  have  been  demon¬ 
strated  recently.5 

Characterization  of  low-frequency  noise  is  needed  to  as¬ 
certain  performance  of  photodetectors.  Noise  measurements 
provide  a  valuable  diagnostic  tool  for  the  evaluation  of  elec¬ 
tronic  and  optoelectronic  devices  and  their  long  term 
performance.6  In  the  case  of  photodetectors,  the  detectivity  is 
limited  by  the  level  of  internal  noise,  usually  expressed  as 
the  noise  equivalent  power  (NEP).  In  this  work,  we  use  low- 
frequency  noise  measurements  of  GaN-based  p-n  junction 
detectors  to  analyze  mechanisms  of  dark  conductivity  and 
device  performance. 

II.  SAMPLE  PREPARATION 

Devices  investigated  here  are  grown  on  sapphire  sub¬ 
strates  by  metalorganic  chemical  vapor  deposition 
(MOCVD).  The  details  of  device  design  and  fabrication  pro¬ 
cedures  are  given  elsewhere.5  The  photodiode  structure  con- 
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sists  of  a  thin  AIN  buffer  layer,  followed  by  a  0.8  /an  thick 
n-type  layer  of  GaN  (doped  cm  3)  and  a  2.2  /zm 

thick  p-type  layer  of  GaN  doped  with  Mg.  The  p -doping 
level  varies  along  the  growth  direction  from  ~  1016  cm"3  at 
the  p-n  junction  interface  to  *^7  X 1017  cm  3  at  the  contact 
surface.  This  results  in  a  p-  'ir-n  structure,  where  tt stands  for 
a  lightly  doped  compensated  p- type  layer  sandwiched  be¬ 
tween  more  heavily  doped  p  and  n  layers.  The  doping  levels 
are  determined  from  capacitance— voltage  (C— V0  and  Hall 
measurements.  Reactive  ion  etching  is  used  to  define  rectan¬ 
gular  mesas  of  individual  diodes  prior  to  Ohmic  contacts 
deposition. 

The  photodiodes  responsivity  reaches  0.14  AAV  at  363 
nm,  and  it  drops  by  more  than  three  orders  of  magnitude 
above  the  cutoff  wavelength.  From  —  2  to  —  20  V  the  pho¬ 
tocurrent  is  almost  independent  of  the  reverse  bias.  The  pho¬ 
todiode  bandwidth  of  ^  32  MHz  is  determined  from  the  mea¬ 
sured  detector  response  to  —  1  ns  long  laser  pulse. 

III.  DARK  CURRENT  MEASUREMENTS 

Reverse  bias  current- voltage  (/-V)  characteristics  of 
GaN  photodiodes  are  measured  in  the  temperature  range 
from  290  to  570  K.  Figure  1  presents  five  of  the  curves, 
taken  on  a  200X  200  /*m2  device.  All  of  the  curves  show  an 
exponential  increase  in  the  dark  current  with  bias  voltage. 
An  increase  in  the  dark  current  with  the  device  temperature, 
for  a  given  value  of  the  reverse  bias,  is  also  approximately 
exponential. 

The  exponential  dependence  of  the  dark  current  on  both 
the  voltage  and  temperature  is  difficult  to  explain  by  conven¬ 
tional  models  for  reverse-bias  conductivity.7  Measurements 
of  devices  with  different  sizes  show  that  the  current  scales 
with  the  device  area,  rather  than  the  perimeter  of  the  active 
region,  thus  ruling  out  any  significant  contribution  of  leakage 
currents  at  the  device  mesa  sidewalls.  Due  to  the  large  band 
gap  of  GaN  and  extremely  small  thermal  excitation  rate,  dif¬ 
fusion  currents  in  the  neutral  region  and  generation- 
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Reverse  bias  (V) 

FIG.  1.  Reverse  bias  branch  of  the  current-voltage  characteristics  for  a 
200X  200  /an2  photodiode  at  five  different  temperatures.  Fitting  results  are 
shown  by  dashed  lines  for  the  low-bias  data  and  by  dotted  lines  for  the 
high-bias  data. 


recombination  (g-r)  currents  in  the  depletion  region  should 
be  negligibly  small.  Even  though  the  room  temperature  dark 
current  for  low  reverse  bias  is  less  than  1  pA,  that  is  still  106 
times  higher  than  the  saturation  current  of  —  2X 10” 18  A  ex¬ 
pected  from  the  thermal  excitation  of  carriers.  The  onset  of 
impact  ionization  in  the  photodiodes  under  test  is  observed 
only  at  a  bias  of  approximately  —42  V,  and  therefore  it  is 
not  expected  to  contribute  significantly  to  the  measured  low 
voltage  /—  V  characteristics.  Finally,  direct  band-to-band  or 
trap  assisted  tunneling  currents  have  exponential  voltage  de¬ 
pendence,  but  should  be  only  weakly  temperature  dependent 

The  conductivity  of  GaN  diodes  can  be  modeled  by  as¬ 
suming  that  the  current  is  due  to  hopping  of  charged  carriers 
via  localized  defect-related  states  (traps)  in  the  depletion  re¬ 
gion.  Hopping  was  studied  extensively  in  amorphous  semi¬ 
conductors,  where  it  was  found  to  be  the  dominant  mecha¬ 
nism  of  current  flow,  especially  at  low  temperatures.8  Hill9 
and  Poliak  and  Riess10  concluded,  using  different  arguments, 
that  for  moderate  electric  fields,  the  current-field  dependence 
is  described  by: 

f  e Fa  lT0\m] 

y=7(0)exp|c^[y]  0) 

where  j  is  the  current  density,  j( 0)  is  the  low-field  current 
density,  T0  is  a  characteristic  temperature  parameter,  T  is  the 
temperature,  k  is  the  Boltzmann  constant,  F  is  the  electric 
field,  e  is  the  electron  charge,  C  is  a  constant  of  the  order  of 
unity,  and  a  is  the  localization  radius  of  the  electron  wave 
function.  At  very  low  fields  hopping  conductivity  cr{F) 
= j(F)/F  is  expected  to  follow  Mott’s  law  for  variable-range 
hopping:11 

<t(0)  a  exp|  —  |  yj l/4] ,  (2) 

where  T{)  is  the  same  parameter  as  in  Eq.  (1).  The  validity  of 
Eq.  (2)  was  experimentally  demonstrated  for  many  amor¬ 
phous  solids.12,13 


We  use  Eq.  (1)  to  fit  the  experimental  I-V  data  for 
reverse  bias  ranging  from  —  1  to  —  6  V.  We  assume  that,  for 
this  bias  range,  the  total  current  through  the  p-n  junction  is 
equal  to  the  saturation  current  originating  in  the  depletion 
region  and  its  density  is  given  by  Eq.  (1).  The  depletion 
region  width  w  was  extracted  from  the  measured 
capacitance-voltage  (C-V)  characteristics.  The  average 
electric  field  dependence  on  the  reverse  bias  V  was  then 
found  from  F=  ( V+  V,)/w,  where  V{  is  the  built-in  junction 
voltage  (estimated  at  —  1  V).  Using  this  relation,  we  convert 
the  measured  current- voltage  data  into  a  current-field  depen¬ 
dence,  at  each  temperature,  which  can  be  fit  directly  with 
Eq.  (1). 

The  best  fit  to  the  experimental  data  for  all  temperatures 
is  shown  in  Fig.  1  by  dashed  lines.  At  each  temperature,  two 
fitting  parameters  are  used:  zero-field  current  density  j( 0), 
and  r=Ca(T0IT)m,  which  has  a  physical  meaning  of  the 
characteristic  jump  length.  The  temperature  dependence  of 
y(0),  fitted  with  Eq.  (2),  is  used  to  extract  the  parameter  T0 
of  1.16X1010  K.  This  value  of  T0  is  approximately  three 
orders  of  magnitude  higher  than  the  107  usually  ob¬ 
served  in  amorphous  semiconductors  such  as  germanium, 
silicon,  or  heavily  implanted  GaAs.14  The  high  value  of  T0 
explains  the  strong  temperature  dependence  of  the  reverse- 
bias  current  in  our  photodiodes.  Assuming  a  localization  ra¬ 
dius  of  a  =  10  A,15  for  a  known  T0,  we  obtain  C=0.4  from 
the  temperature  dependence  of  the  parameter  r,  a  median 
value  between  C=0.8  obtained  by  Hill9  and  C=0.17  ob¬ 
tained  by  Poliak  and  Riess.10  Finally,  the  density  of  states, 
per  unit  energy,  taking  part  in  hopping  conduction  is  esti- 
mated  from  N~W(kT0a3),'6  as  -  1.8X  1016  cm_3eV_l. 

It  is  surprising  to  find  in  a  crystalline  semiconductor  the 
type  of  conductivity  previously  observed  in  amorphous  sol¬ 
ids  or  at  very  low  temperatures.  We  believe  that  the  reason 
for  this  is  high  dislocation  density  in  epitaxial  GaN.  For  a 
typical  dislocation  density  of  108  cm-2,  there  are  on  average 
40  000  dislocations  in  a  200X  200  /m2  device.  The  dangling 
bonds  at  dislocation  boundaries  are  likely  candidates  for  the 
equivalent  of  localized  states  taking  part  in  hopping 
conduction.17 

Equation  (1)  is  expected  to  be  valid  only  for  fields  F 
<^{2kTlea)^5  X 105  V/cm,  at  room  temperature.  It  is  evi¬ 
dent  from  Fig.  1  that  the  slope  of  measured  I-V  character¬ 
istics  changes  noticeably  for  reverse  bias  voltages  higher 
than  ~7  V,  which  corresponds  to  the  electric  field  of  ~2.2 
X  10s  V/cm.  We  believe  that  at  higher  bias  the  diode  con¬ 
ductivity  behavior  can  be  explained  by  a  combination  of 
hopping  and  Poole— Frenkel  effect,  or  field -assisted  thermal 
ionization  of  carriers  from  traps  in  the  depletion  region. 
Similar  behavior  was  previously  observed  in  amorphous 
germanium.8  For  sufficiently  high  electric  fields,  a  single  hop 
can  take  a  charge  carrier  into  the  conduction  band,  instead  of 
taking  it  to  the  neighboring  localized  state.  Poole-Frenkel 
conduction  was  observed  in  thin  films  of  various  wide  band- 
gap  materials,  including  AIN  and  diamond.18,19  It  was  shown 
previously  that  the  field  dependence  of  the  current  density  is 
given  by:20 
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FIG.  2.  Measured  dark  current  noise  spectra  and  the  corresponding  fits  at 
five  different  temperatures. 
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where  j0  is  the  low-field  current  density,  ySpp  is  the  Poole- 
Frenkel  constant,  e0  is  the  vacuum  permittivity,  and  e  is  the 
high-frequency  dielectric  constant  of  the  material.  For 
GaN,  taking  e(°°)  =  5.35,  we  calculate  /Spf=3.3 

xio-4  «vrl(2  cmin. 

The  curves,  calculated  according  to  Eq.  (3),  shown  by 
dotted  lines  in  Fig.  1,  demonstrate  an  excellent  agreement 
with  the  experimental  high  bias  data.  At  all  temperatures,  the 
best  fit  is  obtained  for  /3PF=4.5Xl(T4  eV  V~m  cmm.  It  is 
interesting  to  note  that  for  high  bias,  Eq.  (1)  predicts  currents 
much  higher  than  those  observed  experimentally.  This  can  be 
explained  assuming  that  the  width  of  the  energy  band  for 
localized  states  near  the  Fermi  level  is  limited  to  a  certain 
value  of  A E.  The  enhancement  of  hopping  conduction  with 
electric  field  will  then  saturate  when  the  energy  acquired  by 
a  charge  carrier  moving  against  the  field  eFR ,  where  R  is 
the  length  of  the  jump,  becomes  comparable  to  A E,  At  room 
temperature,  for  F=2.2X  10s  V/cm,  we  estimate  A E 
«0.1  eV.  This  width  is  typical  of  dislocation-related  energy 
bands.17  Consequently,  we  obtain  7V'=WA£**2 
X  1015  cm"3  for  the  spatial  density  of  hopping  conduction 
centers.  Another  conclusion  is  that  even  at  high  bias,  hop¬ 
ping  is  responsible  for  a  significant  part  of  the  diode  current, 
which  can  explain  the  experimentally  estimated  /9PF  being 
~  1 .4  times  its  theoretical  value. 


IV.  NOISE  MEASUREMENTS:  REVERSE  BIAS 

The  noise  characteristics  of  GaN  photodiodes  are  mea¬ 
sured  in  the  frequency  range  of  1  Hz- 100  kHz  using  a  low 
noise  current  preamplifier  and  a  fast  Fourier  transform  (FFT) 
spectrum  analyzer.  The  bias  is  varied  from  “5  to  -30  V, 
and  the  sample  temperature  from  298  to  523  K.  The  equip¬ 
ment  is  calibrated  by  measuring  the  room  temperature  noise 
current  of  a  conventional  600  ft  resistor.  The  noise  floor  of 
the  setup,  determined  by  the  noise  parameters  of  the  preamp¬ 
lifier,  is  ~  10” 25  A2/Hz.  The  experimental  results  for  - 10  V 
bias,  and  different  temperatures,  are  presented  in  Fig.  2.  At 
elevated  temperatures  or  high  values  of  the  reverse  bias. 


FIG.  3.  Bias  and  temperature  dependence  of  the  proportionality  coefficient 
relating  the  dark  current  and  noise  power. 


where  the  dark  current  of  the  photodiode  increases  above 
~  1  nA,  the  appearance  of  Ilf  noise  is  observed.  All  of  the 
measured  noise  spectra  satisfy  the  usual  relationship: 


(4) 


where  Sn  is  the  spectral  density  of  the  noise  current.  Id  is  the 
total  current,  /  is  the  frequency,  and  Jq  T  are  fitting 
parameters.  The  value  of  y  is  found  to  vary  from  1.0  to  1.1. 
The  values  of  Sq  are  plotted  in  Fig.  3  as  a  function  of  bias 
and  temperature.  The  parameter  s0  decreases  slightly  with 
increasing  bias  and  rapidly  with  increasing  temperature. 

Two  models  of  the  1  If  noise  in  hopping  conductivity 
have  been  proposed,  both  of  which  are  related  to  the  pres¬ 
ence  of  localized  states,  separated  from  each  other  by  dis¬ 
tances  longer  than  the  typical  electron  jump  length.  Such 
sites  are  not  incorporated  into  the  so-called  critical  current- 
carrying  network,  but  function  as  “dead  ends,  slowly  ex¬ 
changing  electrons  with  sites  within  the  network.  In  the  num¬ 
ber  fluctuation  theory,21  this  leads  to  slow  fluctuations  in  the 
total  number  of  electrons  taking  part  in  hopping  conduction. 
The  “mobility”  fluctuation  theory22  assumes  that  an  electron 
trapped  at  a  “dead  end”  changes  the  distribution  of  charge 
in  the  surrounding  material  and  shifts  the  energy  levels  of 
neighboring  traps.  This  leads  to  slow  fluctuations  in  the  ef¬ 
fective  “mobility”  of  hopping  electrons. 

The  “mobility”  fluctuation  theory  predicts  a  very  weak 
temperature  dependence  of  the  low-frequency  hopping  noise 
and  it  does  not  appear  applicable  in  our  case.  The  number 
fluctuation  theory  results  in  the  frequency  dependence  of  the 
noise  spectral  density  given  by  a  well-known  Hooge 
formula:23 


(5) 


where  N  is  the  average  number  of  electrons  taking  part  in  the 
conduction  process  and  the  predicted  value  of  the  Hooge 
parameter  a  is  of  the  order  of  unity.  We  can  roughly  estimate 
the  value  of  a  for  our  photodiodes  from  the  following  con- 
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FIG.  4.  Normalized  noise  spectra  at  2  V  forward  bias  and  different 
temperatures. 


FIG.  5.  Arrhenius  plot  for  the  trap  center  lifetime  obtained  from  the  mea¬ 
sured  noise  spectra. 


sideration.  From  Eqs.  (4)  and  (5)  we  have  oc=SqN.  At  room 
temperature  and  — 10  V  reverse  bias,  the  measured  value  of 
the  parameter  sQ  is  ~1.3X1(T5.  At  the  same  temperature 
and  bias  the  dark  current  of  the  photodiode  is  ~  1  nA,  or, 
equivalently,  the  sample  resistance  is  J^IXIO*0  ft.  The 
number  of  electrons  N  can  be  found  from  R  —  L?/(e/iN)t 
where  the  length  L  is  taken  to  be  equal  to  the  depletion  layer 
width  w,  obtained  from  C-  V  measurements.  The  effective 
mobility  /z  is  estimated  from  the  theory  of  hopping 
conduction9  using  parameters  obtained  from  the  fit  to  the 
/-V  data  at  room  temperature.  We  obtain  yi— 5 
X  10~6  cm2  V"1  s~\  JV~2.3X  105,  and  finally,  a~3,  con¬ 
sistent  with  the  predictions  of  the  number  fluctuation  theory. 
Both  the  bias  and  temperature  dependence  of  the  parameter 
sq  can  be  explained  by  the  change  in  the  number  of  hopping 
electrons  N.  With  increasing  bias,  N  goes  up  due  to  increas¬ 
ing  width  of  the  depletion  layer,  and  with  temperature,  due  to 
decreasing  length  of  a  typical  jump. 

V.  NOISE  MEASUREMENTS:  FORWARD  BIAS 

Under  a  forward  bias,  the  built-in  electric  field  de¬ 
creases,  and  hopping  is  no  longer  a  dominant  conductivity 
mechanism.  The  photodiode  current  is  dominated  by  diffu¬ 
sion  currents  in  the  neutral  region.  As  a  consequence,  we  do 
not  observe  the  1  If  current  noise.  Instead,  most  devices  show 
a  very  clear  g-r  noise  component  with  a  Lorentzian  spec¬ 
trum  of  the  form: 

5"(<W)=A  1  +(<UTo)i’  (6) 

where  A  is  a  constant  and  r0  is  the  characteristic  time  of  the 
generation-recombination  process.  For  the  accurate  determi¬ 
nation  of  r0  from  experimental  data,  it  is  convenient  to  mul¬ 
tiply  the  measured  noise  density  S„/Ij  by  (o—2irf  and  plot 
it  on  a  linear  scale  as  a  function  of  frequency  f.u  The 
Lorentzian  then  appears  as  a  symmetric  peak  at/0,  and  To  is 
found  from  t0=1/2it/o. 

In  Fig.  4,  we  plot  the  frequency  dependence  of  the  nor¬ 
malized  noise  density  at  three  different  temperatures,  mea¬ 


sured  on  a  200X200  /xm2  device  at  a  forward  bias  of  +2  V. 
It  is  clear  from  Fig.  4  that  /o  is  increasing  and,  therefore.  To 
is  decreasing  with  temperature.  For  a  thermally  activated  g-r 
process  the  characteristic  time  r0  is  expected  to  follow  an 
Arrhenius  equation: 


7o=*bo  exp 


(7) 


where  Too  *s  a  constant  and  €  is  die  activation  energy.  The 
value  of  e  can  be  estimated  from  the  slope  of  log(r(172)  plot¬ 
ted  against  the  inverse  temperature.25  The  Arrhenius  plot  for 
our  device  is  presented  in  Fig.  5.  We  obtain  e— <0.49  eV,  in 
agreement  with  the  expected  activation  energy  for  centers 
associated  with  Ga  antisite  defects  in  GaN.26 


VI.  NOISE  MEASUREMENTS  UNDER  ILLUMINATION 

We  have  also  performed  noise  measurements  on  GaN 
photodiodes  illuminated  with  a  Xe  arc  lamp.  The  photodiode 
bias  was  kept  at  -5  V,  and  the  light  intensity  was  adjusted 
by  changing  the  lamp  current.  Unfortunately,  at  low  frequen¬ 
cies  the  measured  noise  was  dominated  by  the  noise  of  the 
lamp  itself,  so  we  do  not  have  experimental  evidence  for  the 
presence  or  absence  of  the  1  If  noise  in  the  photocurrent. 
Above  — 10  kHz,  the  measured  noise  spectrum  was  com¬ 
pletely  flat  The  corresponding  noise  current  density  is  plot¬ 
ted  in  Fig.  6  as  a  function  of  photocurrent  lp .  Figure  6  shows 
that  the  experimental  points  lie  very  close  to  the  line  of 
yj2elp,  expected  from  the  shot  noise  theory  27  We  conclude 
that  at  least  at  these  frequencies,  no  other  significant  noise 
sources  are  present  in  photoconduction. 

VII.  DEVICE  PERFORMANCE 

At  room  temperature  and  for  reverse  bias  values  lower 
than  10  V,  the  dark  current  noise  spectra  disappear  below  the 
noise  floor  of  our  measurement  setup.  However,  the  device 
performance  can  be  estimated  by  extrapolating  noise  data 
taken  at  higher  biases.  The  dark  current  at  -  3  V  is  2.7  pA, 
resulting  in  Sn  (at  1  Hz)  =7.3X10  29  A2/Hz.  At  this  bias 
the  magnitude  of  the  1  If  noise  for  frequencies  above 


r 


2146  J.  Appl.  Phys.,  Vol.  83,  No.  4, 15  February  1998 


Kuksenkov  et  a!. 


—  85  Hz  will  be  lower  than  that  of  the  shot  noise, 
S^MXHT31  A2/Hz.  The  corresponding  NEP  (/ 
^  100  Hz)  is  6.6X10”15  W/Hz1/2,  comparable  to  the  best 
data  reported  for  ultraviolet-enhanced  silicon  photodiodes. 

Our  analysis  shows  that  leakage  current  and  noise  pa¬ 
rameters  of  GaN  photodetectors  are,  at  present,  limited  by 
the  quality  of  epitaxial  material.  Reducing  the  density  of  dis¬ 
locations  and  point  defects  would  eliminate  hopping,  and 
significantly  reduce  the  dark  current  and  noise  intensity.  In 
addition,  since  hopping  conduction  is  greatly  enhanced  by 
the  electric  field,  we  expect  that  the  device  performance  will 
be  improved  by  fabricating  p-i-n  diodes  with  a  wider  intrin¬ 
sic  region. 

VIII.  CONCLUSION 

We  have  studied  the  origins  of  dark  conductivity  and 
low-frequency  noise  in  the  dark  current  of  GaN  p-n  junction 
photodetectors.  Hopping  via  defect-related  states  in  the 
depletion  region,  thought  to  be  related  to  dislocations,  is 
found  to  be  the  dominant  current  mechanism.  The  density  of 
localized  states  in  the  forbidden  gap  is  estimated  at  2 
X  1015  cm-3.  Different  noise  sources  are  identified  in  the 
dark  current  under  forward  and  reverse  bias,  and  under  illu¬ 
mination.  Under  reverse  bias,  the  dark  current  noise  is  of  the 
Mf  character  and  obeys  the  Hooge  relation  with  0^3.  Under 
forward  bias,  generation-recombination  noise  related  to  the 
trap  level  with  the  activation  energy  of  0.49  eV  is  observed. 


Under  illumination,  the  photocurrent  noise  is  found  to  be 
dominated  by  the  shot  noise.  The  noise  equivalent  power  of 
a  200X  200 /an2  device  at  -3  V  bias  is  estimated  at  6.6 
XKT15  W/Hzl/2. 
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